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Abstract 


We  summarize  progress  in  our  analysis  of  proposed  communication  schemes  for  hypersonic  ve¬ 
hicles  in  flight.  Investigations  into  the  electron-acoustic  wave  (EAW)  communication  scheme  have 
yielded  a  dispersion  analysis  indicating  that  EAW  modes  can  be  generated  and  propagate  in  the 
plasma  layer.  Furthermore,  we  have  demonstrated  that  these  modes  can  be  coupled  to  electro 
magnetic  waves  at  the  plasma  layer  boundarv  for  the  assumption  of  a  sharp  density  boundary.  A 
detailed  analysis  of  the  wave  transformation  coefficients  at  this  boundary  layer  is  presented.  We 
also  discuss  the  results  of  ID  particle-in-cell  simulations  of  a  plasma- filled  cross-field  diode.  These 
simulations  are  intended  to  model  some  of  the  critical  physical  processes  of  the  ReComm  scheme 
for  communications  through  the  plasma  sheath  surrounding  a  hypersonic  vehicle  during  re-entry. 
The  ReComm  scheme  utilizes  externally  generated  crossed  electric  and  magnetic  fields  to  induce 
charged-particle  drifts.  This  spatially  redistributes  the  plasma  and  enables  propagation  of  electro¬ 
magnetic  waves  through  the  plasma  layer.  We  demonstrate  that  the  time  required  to  significantly 
alter  the  plasma  density  profile  in  the  diode  corresponds  roughly  to  the  ambipolar  diffusion  time. 
When  a  voltage  is  applied,  the  plasma  center-of-mass  shifts  from  the  center  of  the  diode.  But 
whether  that  shift  is  towards  the  cathode  or  anode  depends  sensitively  on  the  relative  mobility 
of  the  plasma  ions  and  electrons.  Finally,  we  summarize  the  derivation  of  an  analytic  model  of 
whistler  wave  excitation  by  an  external  dipole  antenna.  In  addition,  a  dispersion  relation  is  derived 
to  estimate  the  plasma  heating  in  the  sheath  due  to  plasma  waves  excited  by  the  antenna. 
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I.  INTRODUCTION 


Hypersonic  vehicles  traveling  in  the  upper  atmosphere  at  speeds  greater  than  Mach  10 
(7,000  mph  at  33  km  altitude)  generate  plasma  that  disrupt  or  prevent  communications  over 
conventional  radio-frequency  channels.  This  phenomenon  manifests  itself  most  famously  as 
the  well-known  communications  “blackout’’  period  during  space  vehicle  reentry  into  the 
atmosphere  [1-3],  and  has  been  studied  (off  and  on)  since  the  Gemini  and  Apollo  space 
programs.  For  sustained  hypersonic  flight  in  the  atmosphere,  the  communications  blackout 
will  persist  for  almost  the  entire  flight. 

The  plasma  boundary  layer  is  formed  around  a  hypersonic  vehicle  in  the  bow  shock  and 
to  some  extent  by  the  collision  of  neutral  gas  particles  with  material  ablated  from  the  surface 
of  the  hypersonic  vehicle  itself.  In  the  shock,  the  air  is  heated  to  such  a  high  temperature 
that  ionizing  collisions  between  neutral  particles  occur.  Estimates  of  electron  density.  np 
(~  1012  cm-3),  temperature  (~  0.5  eV)  and  electron-neutral  collision  frequencies  (2.5  x  109 
42.5  x  1011  s  1 )  have  been  obtained  for  high-Mach-number  reentry  vehicles  at  altitudes 
between  100,000-150,000  ft  (30-45  km)  [4,  5].  The  thickness  of  the  plasma  layer  (defined  as 
the  region  where  the  electron  density  exceeds  ~  109  cm-3)  associated  with  these  parameters 
is  on  the  order  of  6  cm,  with  90%  of  the  peak  density  out  to  a  width  of  about  2  cm.  The  free 
electrons  in  the  plasma  layer  attenuate  radio-frequency  (RF)  waves  both  through  reflection 
and  resistive  absorption,  and  generate  RF  noise.  At  sufficiently  large  antenna  power  levels, 
the  RF  itself  can  cause  further  ionization  of  the  air  to  occur. 

Acceptably  low  attenuation  of  electromagnetic  (EM)  wave  propagation  in  a  collisionless 
plasma  is  limited  to  frequencies  above  the  cutoff  value  fc  =  u:pfJ‘2z.  where  cope  is  the  electron 
plasma  frequency  (see,  for  example.  Ref.  [6]).  For  np  -  1012  cm-3,  the  fc  ~  10  GHz,  which 
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is  at  the  high  end  of  the  radio  wave  frequency  spectrum  and  well  into  the  conventional 
microwave  portion  of  the  EM  spectrum.  In  addition,  this  plasma  cutoff  frequency  is  above 
the  spectrum  allocated  (in  some  cases  by  international  treaty)  for  flight  test  and  evaluation 
telemetry  and  above  frequencies  supported  by  existing  infrastructure.  We  note  that  Global 
Positioning  System  signals  are  broadcast  in  the  L-band  region  of  the  spectrum,  and  con¬ 
sequently  are  also  heavily  attenuated  by  the  plasma  sheath.  As  a  result,  communication 
approaches  within  the  standard  telemetry  spectrum  are  inadequate  for  hypersonic  vehicles 
using  conventional  techniques  (see,  for  example.  Refs.  [6,  7]). 

A  theoretical  research  program  is  underway  to  assess  concepts  for  real-time  telemetric 
communications  through  plasma  boundary  layers  surrounding  hypersonic  vehicles  in  flight. 
This  research  includes  detailed  analytic  analysis  and  computational  modeling  of  several 
proposed  communication  schemes.  In  our  previous  report  [8],  we  summarized  our  research 
progress  on  three  such  schemes:  a  nonlinear  three- wave  interaction  technique  [9-11],  an 
electron-acoustic  wave  (EAW)  propagation  scheme,  and  a  magnetic  window  scheme  for  EM 
wave  propagation  through  a  collisional  plasma  slab.  Here  we  summarize  recent  progress 
in  our  analysis  of  the  EAW  communication  scheme  (Sec.  II).  In  a  separate  report  [12], 
we  present  updated  2D  numerical  simulations  of  the  magnetic  window  scheme  and  new’ 
3D  simulation  results.  These  new’  simulations  demonstrate  the  ability  to  model  detailed, 
complex  antenna  structures  which  is  an  important  aspect  to  all  the  communication  schemes 
presently  under  analysis. 

In  Sec.  Ill,  we  present  the  results  of  ID  particle-in-cell  simulations  of  a  plasma-filled  cross- 
held  diode.  These  simulations  are  intended  to  model  some  of  the  critical  physical  processes 
of  the  ReComm  scheme  for  communications  through  the  plasma  sheath  surrounding  a  hy¬ 
personic  vehicle  during  re-entry.  The  ReComm  scheme  utilizes  externally  generated  crossed 
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electric  and  magnetic  fields  to  induce  charged-particle  drifts.  This  spatially  redistributes 
the  plasma  and  enables  propagation  of  electromagnetic  waves  through  the  plasma  layer.  We 
demonstrate  that  the  time  required  to  significantly  alter  the  plasma  density  profile  in  the 
diode  corresponds  roughly  to  the  ambipolar  diffusion  time.  When  a  voltage  is  applied,  the 
plasma  center-of-mass  shifts  from  the  center  of  the  diode.  But  whether  that  shift  is  towards 
the  cathode  or  anode  depends  sensitively  on  the  relative  mobility  of  the  plasma  ions  and 
electrons. 

In  Sec.  IV,  we  present  analytic  modeling  of  the  whistler  wave  communication  scheme. 
The  impact  of  the  antenna  radiation  field  on  the  plasma  sheath  properties  is  assessed.  A 
summary  of  our  progress  on  the  various  communication  schemes  is  given  in  Sec.  V. 

II.  ELECTRON  ACOUSTIC  WAVE  TRANSFORMATION 

One  scheme  for  communication  through  a  plasma  sheath  is  based  on  the  possibility  of 
propagating  EAWs  with  frequencies  well  below  the  local  plasma  frequency  (~  9  GHz)  in  a 
two-electron-temperature  plasma.  A  small  hot  electron  population  is  produced  by  injection 
of  an  energetic  electron  beam.  Over  a  certain  range  of  wavelengths,  the  hot-electron  com¬ 
ponent  of  the  plasma  can  short  out  the  electric  fields  that  produce  plasma  oscillations  in 
the  cold  plasma,  thereby  reducing  the  frequency  of  these  oscillations  below  the  cold  plasma 
frequency.  It  may  be  possible  to  excite  these  electrostatic  oscillations  with  an  antenna  and 
have  them  couple  to  electromagnetic  RF  oscillations  outside  the  sheath.  Related  research 
was  carried  out  in  the  CHARGE  2B  ionospheric  rocket  experiments  [13-15].  For  sustained 
hypersonic  flight,  a  comprehensive  analysis  is  required  to  estimate  the  feasibility  of  this 
concept. 
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A.  EAW  to  EM  Wave  Transformation  at  the  Sheath  Boundary 


We  use  the  hydrodynamic  approximation  to  describe  the  transformation  of  an  EAW  into 
an  EM  wave  at  the  plasma  sheath  boundary.  This  transformation  is  essentially  the  conver¬ 
sion  of  longitudinal  EAWs  propagating  inside  the  plasma  layer  into  EM  waves  propagating 
in  the  atmosphere. 

We  retain  the  two  temperature  plasma  description,  with  the  cold  and  hot  electron  tem¬ 
peratures  denoted  as  Tc  and  Th,  respectively.  The  total  electron  plasma  density  is  the  sum 
of  the  hot  and  cold  electron  components,  Ne  =  nc  +  rih .  Density  perturbations  under  the 
influence  of  the  wave  fields  are  denoted  as  nc  =  noc  +  finc  and  ///,  =  ni)h  +  dnh,  and  electron 
velocity  perturbations  are  denoted  as  6vc  and  <5tv  We  consider  the  transformation  of  an 
incident  EAW  with  frequency 


(1) 


propagating  in  the  plasma  layer  into  an  EM  wave  propagating  in  air, 


(2) 


[The  notation  used  here  follows  Ref.  [8];  e.g.,  ujvc  is  the  plasma  frequency  of  the  cold  electrons 
and  rph  is  the  Debye  radius  of  the  hot  electrons.] 

We  assume  that  in  the  case  of  EAWs  propagating  in  the  plasma  sheath,  in  the  vicinity 
of  the  boundary  there  are  two  types  of  waves;  EM  waves  with  amplitude  exponentially 
decreasing  with  increasing  plasma  density,  and  electrostatic  EAWs  propagating  towards  the 
boundary.  In  addition,  we  also  take  into  account  reflected  EAWs  and  EM  waves,  propagating 
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FIG.  1  Incident  longitudinal  and  reflected  EAWs  (in  blue  with  subscript  s)  and  the  transformed 
and  reflected  EM  waves  (in  red  with  subscript  t)  of  P-polarization  ( Ex ,  Ez,  Hy). 

in  the  air  (see  Fig.  1). 

In  order  to  derive  equations  for  the  transformation  of  EAWs  at  the  boundary  to  P- 
polarized  EM  waves,  we  begin  with  Maxwell’s  equations: 


V  x  E  =  — 


ldH 
c  dt 


(3) 


(4) 


divE  =  — 47re  (nlc  +  n\h)  ■ 


The  linearized  equation  of  motion  for  the  cold  electrons  can  be  written  as 


(6) 
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where  vec  is  the  cold  electron-neutral  collision  frequency,  and  the  (first-order)  cold  electron 
fluid  pressure  is  expressed  as 

VPic  =  70  mV£.Vnlc, 


with  70  representing  the  adiabatic  constant,  and 


=  s/%R. 


(Subscripts  prefixed  by  “0”  and  “1”  denote  the  usual  zeroth  and  first  order  terms,  respec¬ 
tively,  in  the  linear  expansion.)  The  linearized  continuity  equations  for  the  cold  electrons 
has  the  form 

+  n0rdmJlc  =  0.  (7) 

Similarly  for  the  hot  electrons,  the  linearized  equations  of  motion  and  continuity  are 


dvi  h 
df 


e  -  l  _ 

P  ^  Pi  h 

m  mriQC 


(8) 


dnlh 

dt 


+  noh'V  •  Vih  —  0. 


(9) 


with 


VPi*  =  7 0mV£hVnlh, 


and 


VTh  =  VThJ™- 


Equations  (4)  -  (9)  can  be  used  to  describe  both  (a)  the  transformation  of  an  incident 
EM  wave  into  an  EAW  inside  the  plasma  sheath  as  well  as  (b)  the  transformation  of  an 
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EAW  incident  on  the  plasma/air  boundary  into  an  EM  wave.  Here  we  focus  on  these 
transformations  at  a  relatively  sharp  boundary  of  a  plasma  sheath.  If  L  is  the  characteristic 
distance  over  which  the  plasma  density  is  decreasing  from  some  finite  value  to  zero  and  5  is 
the  skin  depth  for  EM  held  penetration  into  the  plasma,  then  the  condition  for  which  the 
plasma  boundary  can  be  defined  as  sharp  is  simply  S  L. 

Applying  the  curl  operator  to  Eq.  (4)  and  using  Eqs.  (3)  and  (7),  the  following  equation 
for  the  y-component  of  magnetic  held  in  plasma  can  be  obtained: 

0,  (10) 


where  ktz  is  the  transformed  EM  wave  vector  (note  that  ktz  —  kSz)i  penetrating  into  the 
plasma  sheath  and 


(11) 


Taking  into  account  that  ec  <  0,  we  see  from  Eq.  (11)  that  the  magnetic  held  associated 
with  the  penetrated  EM  wave  will  exponentially  decrease  as  the  wave  propagates  in  the  x 
direction. 

Now  wre  can  obtain  expressions  for  the  perturbed  electron  velocities  V\c  and  Vy,  for  use 
in  Eq.  (4).  From  Eqs.  (6)  and  (7),  it  is  straightforward  to  obtain  the  following  expression: 


1  + 


7»V?C 


w  (uj  +  iVec) 


Vdiv 


Vic  =  -t 


miv  (uo  +  iuec) 


(12) 


Taking  into  account  that  for  the  cold  plasma  electrons: 

E »  .7oV^c  VdivE, 

U![U>  +  IV, eC) 


(13) 
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we  can  rewrite  Eq.  (12)  for  V\c  as 


i— _ 1 _ B  +  i- . 

m  u  (u  +  ivec)  mu{ u>  +  ivec) 


VdivE. 


Introducing  the  linear  operator,  L/,, 


(14) 


Lh  —  1  + 


7o  Vti 


Th 


uj{uj  +  iveh) 


A, 


(15) 


and  noting  that  for  the  hot  electrons  the  second  term  in  (15)  is  larger  than  the  first  one, 
from  Eqs.  (8)  and  (9)  we  can  obtain  the  following  expression  for  the  linearized  speed  of  the 
hot  electrons  in  the  plasma  sheath: 


1 


+  iueh) 


1  + 


7oV£ 


Th 


U>(uJ  +  IVeh) 


A 


E  +  i- 


701?, 


Th 


mu(uj  +  wehy 


:VV  •  E. 


(16) 


Now,  applying  the  operator  L/,  to  Eq.  (4)  and  using  Eqs.  (14)  and  (16),  we  obtain  the 
following  equation: 


V  x 


( Lhti )  = 


.UJ 


UJ, 


ph 


+  iveh) 


oj: 


LhE  -  i- 


ph 


loVl 


C  Uj{uj  +  iueh) 


— — ^VdivE.  (17) 


When  applied  to  the  transverse  EM  held,  L^  —  1.  Rewriting  Eq.  (3)  as 


i J  dEy 

i—Hy  =  ikzEx  — jA 
c  ax 


(18) 
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Eq.  (17)  can  be  expressed  as: 


7oV£ 


d2Ex  u 

- -  H — 

'c((jJ  +  IVeh)  dx 2  C 


Th 


ph 


UJ  {iO  +  ~ 


(l-  l°Vlhk‘‘  ,)  E, 
\  U  (u  +  iveh)  ) 


(19) 


=  1- 


no^Th^lh  |  A  H 

C?(u>  +  iveh)2  V ’ 


a; 


ph 


jjJ  (u)  +  iveh) 


—  er 


7qV&  9^  t  _ 

(w  +  il/eh)  9a:2  [(w  +  il/efc) 


+ 


ph 


-er  \  lj 


(w  +  ii/gfc)/ 


£2 


(20) 


9Er  _  ^ .  9//,y 
?  <2a;(a;  +  ii4,h)2  9a:  ?r  9a:  ' 


The  system  of  equations  (10),  (19)  and  (20)  describes  the  transformation  of  an  EAW  incident 
onto  the  plasma/air  interface  from  within  the  plasma  into  an  EM  wave.  In  addition,  this 
same  system  of  equations  describes  the  transformation  of  EM  wave  incident  from  the  air 
into  EAW  propagating  inside  the  plasma  sheath. 


B.  Solution  Set 


Transformation  of  EAW  waves  on  the  sheath  boundary  can  be  described  by  the  equations 
(10),  (19),  and  (20).  Solution  of  these  equations  can  be  presented  in  the  form: 


Hlv  =  — £—etoiGEx  exp  {-k  ltxx) ,  (21) 

CK\zt 

EXx  =  Ex  exp  (—iqxx)  -  GEX  exp  (— «# Uxx) ,  (22) 

£i2  =  -Ei—  exp(-iqxx)  +  iG^-Ex  exp  (-n'xtxx) ,  (23) 

<7i  kxtz 

to 

H2y  =  — —etotGE2 exp(—n  xtxx),  (24) 

CKxzt 
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FIG.  2:  Incident  longitudinal  and  reflected  EAW  waves  (in  blue  with  index  s)  and  transformed 
and  reflected  electromagnetic  waves  (in  red  with  index  t)  of  P-polarization  (Ex,  Ez,  Hy). 


E2x  =  E2exp(iqxx)  -  GE2ex p{-k'Uxx), 


^l.s; 


ltx 


E\z  =  E2—  exp(iqix)  +iG- — E2exp(-K'Uxx). 

k\tz 


In  the  solution  set  (21)  -  (26),  the  following  notations  are  used  (see  Fig.  2): 


etot  —  1  ~ 


cu: 


pc 


u, 


ph 


lo(uj  +  ivec)  u:(uj  4-  iveh) 


(25) 

(26) 


(27) 


the  constant  G  is  defined  as: 


E\;2tx  —  GE  i,2sx- 


where  for  El  2  we  have 


E  1,2  —  \EU2sx\, 


,  _  w 
K  ltx  —  — 


C 


+  sin2  03, 


(28) 


(29) 

(30) 
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where  63  is  the  angle  of  the  transmitted  EM  wave  with  respect  to  the  plasma-vacuum 
interface  normal,  and  we  make  use  of  the  conditions 


U) 

Msz  Mix  ^2 sz  ^2 tz  k-3tz  sin  #3. 

C 


(31) 


Finally  for  q\  we  have 

Qi  =  l^’isxl-  (32) 

To  find  G  we  use  the  boundary  condition  for  the  speed  of  cold  and  hot  electrons;  that  tile 
normal  component  of  hydrodynamic  speed  of  hot  and  cold  electrons  is  zero  on  the  boundary: 


Vhx(O)  =  0,  Ucx(0)  =  0. 


(33) 


This  means  that  the  electrons  are  elastically  reflected  from  the  boundary.  Using  the  equa¬ 
tions  of  motion  for  the  cold  and  hot  electrons  and  the  Poisson  equation,  we  write  the  following 
condition  for  the  normal  component  of  the  electric  held  on  the  plasma  sheath  boundary: 


at  x  =  0, 


(34) 


where  for  Rp  is 


i  _  i  1 

p>2  „2  «2 

nn  '  Dc  '  Dh 

Boundary  condition  (31)  allowrs  us  to  determine  the  constant  G,  which  connects  E\xt  with 


E\sxi 


G=  1  + 


m 


D 


Dh 


(35) 
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The  boundary  conditions  for  the  magnetic  field  are 


H  i„(0)  =  H7v(  0)  =  H„{  0). 


(36) 


and  the  electric  held  are 


Eu(0)  +  E2z(0)  =  E3z(0). 


(37) 


Using  Eqs.  (21)  -  (26),  we  obtain  the  following  expressions  for  the  magnitudes  of  reflected 
EAW  ( £2 )  and  transmitted  EM  waves  (£'3): 


Z\  —  cos  63 
cos  03  -  Z2 ' 


(38) 


£3  = 


From  the  boundary  conditions,  the  impedances  Z\  and  Z2  can  be  obtained: 


(39) 


£iz(0)  _  .ck'ux  u sin2  O3 

Hiy{  0)  ujetot  cqiGetot 


(40) 


£^(0)  _  .CK'2tx  u sin2  O3 

H2y(0)  tot  eqiGetot 


(41) 


Now  the  coefficient  of  transformation  in  amplitude  of  an  electron  acoustic  wave  inside  a 
plasma  sheath  into  an  electromagnetic  wave  propagating  in  the  air  can  be  written  as: 


f-1 

1  Vtot\  sin  63 

£1  1 

\cqj 

[u>/  (c^jG)]  sin2  O3  -  |eto*|cos03  +i\J\ec\ 

+  sin2  O3 

In  the  process  of  derivation  of  (42),  we  also  assumed  that  vec  <SC  u). 


(42) 
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The  power  associated  with  the  transformation  of  an  EAW  into  an  EM  wave  can  be  defined 
as  the  ratio  of  the  normal  component  of  the  energy  density  flux  of  an  EM  wave: 

E 2 

Sem  =  c-1  cos  03  (43) 

47T 

to  the  normal  component  of  the  energy  density  flux  of  an  EAW: 


Sleaw  —  Seaw  cos  6\  ~  Seaw , 


(44) 


and  for  Seaw  we  have 


c  „  5  ^ 

SEAW  =  \gJ-^s)  —  . 


(45) 


where  Vg  is  the  wave  group  velocity.  Here  Es  is  the  electrostatic  field  magnitude  for  the 
EAW  and 

„  1  U BC 

f  -  1-1-  -  -  - — 

1  ^  Z,2r2  <2  • 

KsrDh  U 


From  (45)  we  obtain 


*>EAW  — 


1  Upc  El 
ksf Dh  47T 


(46) 


Using  Eqs.  (44)  and  (46),  the  transformation  coefficient  Wy- 


WT  =  167T 


Vrh  Upc  \etot\2  1 


COS  03  sill2  03 


C  Uph  |ec|2  q\r2Dh  [w/  (cqiG) sin2  d3  -  |etot|  cos#3]2  +  (|ec|  +  sin2/?3) 


(47) 


To  obtain  (47),  we  also  used: 


lCcl  ;,2r2 

KsrDh 


1  _ 

^  ' 


(48) 
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Taking  into  account  that  cos#is  ~  1  and  \etot\  ~  |ec|,  we  can  replace  (47)  by: 


Wt  = 


167T  |fc| 


V77 i^pc _ cos  O3  sin2  O3 _ 

c  ^ph  [a ;/  (cqiG) sin2  93  -  \ftot\ cos 03]2  +  (|ec|  +  sin2 03) 


(49) 


In  the  next  section  we  investigate  the  dependence  of  the  transformation  coefficient  Wf 
from  the  propagation  angle  03  of  an  outgoing  electromagnetic  wave  propagating  outside  the 
plasma  sheath. 


C.  Numerical  Evaluation  of  the  Transmitted  Power 


In  the  previous  section,  the  angular  dependence  of  the  transformation  coefficient  Wt 
of  an  EAW  propagating  inside  a  plasma  sheath  into  an  EM  wave  propagating  outside  the 
plasma  sheath  was  obtained  [Eq.  (47)].  Taking  into  account  that  cos 9\s  ~  1  and  |^tot |  ~  |ec|, 
Eq.  (49)  was  obtained.  In  our  analysis  of  this  expression,  we  use  the  following  notations: 

V77 1  Upc  U)  Vih  Upc 

Ui  =  - .  <Z2  —  - •  fl3  —  -  — - =  dlO-2' 

C  L^ph  C  ^ ph 


Ql  —  l^’lsxl  —  kis  COS  6>is  h\3  — 


r“ 

Considering  vec  <?C  u>  and  ueh  <C  u  for  | ,  we  can  simply  wTrite: 


M  — 


Upc  +  Uph 
UJ2 


-  1. 


Now  we  can  rewrite  the  expression  for  the  transformation  coefficient  [Eq.  (49)]in  the  form: 


Wt  =  167r|ec|  0102 


cos  #3  sin2  03 


(a3G~x  sin2  03  -  |et0,|  cos03)  +  (|ec|  +  sin2  03) 


(50) 
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Coefficient  of  transformation 
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FIG.  3:  Dependence  of  the  transformation  coefficient  from  the  propagation  angle  (angle  between 
the  direction  of  the  EM  wave  propagation  and  the  normal  to  a  plasma  sheath)  of  an  outgoing  EM 
wave. 

Equation  (50)  is  numerically  evaluated  (see  Appendix  A)  to  plot  the  dependence  of  Wt  on 
the  EM  wave  transmission  angle  0:i.  We  find  in  Fig.  3  that  the  transformation  coefficient 
is  large  for  large  propagation  angles  of  the  outgoing  EM  wave.  The  transmitted  power  is  a 
maximum  for  propagation  angles  close  to  7t/2  (EM  waves  propagating  almost  parallel  to  the 
plasma  sheath  boundary).  This  is  different  from  the  case  of  Langmuir  wave1  transformation 
into  an  EM  wave  at  the  plasma  boundary,  where  the  transformation  coefficient  is  large  for 
small  propagation  angles  of  the  outgoing  EM  wave.  We  note  that  the  magnitudes  of  these 
transformation  coefficients  are  comparable  in  both  the  Langmuir  wave  and  EAW  cases. 

III.  ID  ANALYSIS  OF  THE  RECOMM  SCHEME 

We  present  the  results  of  ID  particle-in-cell  (PIC)  simulations  of  a  plasma-filled  cross- 
held  diode  using  the  Lsp  code[16].  These  simulations  are  intended  to  model  the  ReComin 
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scheme[17]  for  communications  through  the  plasma  sheath  surrounding  a  hypersonic  vehicle 
during  re-entry.  As  is  well  known,  an  electromagnetic  wave  at  angular  frequency  uj  can 
only  propagate  through  a  thick  collisionless  plasma  with  a  density  low  enough  so  that  the 
plasma  frequency  u>p  oc  \fn  does  not  exceed  u.  In  certain  altitude  ranges  the  density  in  the 
shock  induced  plasma  layer  surrounding  hypersonic  vehicles  is  high  enough  to  attenuate  RF 
signals  in  the  L  band  (/  ~  1  —  3  GHz)  which  is  used  for  telemetry  and  evaluation  (T&E). 

A  simplified  picture  of  the  ReComm  scheme  can  be  seen  in  Fig.  4.  The  antenna  aperture 
is  placed  in  between  a  pair  of  biased  electrodes.  An  electromagnet  provides  a  magnetic 
held  normal  to  the  antenna  aperture  throughout  the  region  between  the  electrodes.  We 
consider  in  this  report  the  “plasma-optic”  regime  described  by  Keidar  et  nl.  [18]  in  which 
electrons  are  magnetized  and  ions  are  unmagnetized.  The  separation  distance  between  the 
plates,  L,  is  assumed  to  be  large  compared  to  the  electron  Larmor  radius,  but  small  with 
respect  to  the  ion  Larmor  radius.  In  this  case,  if  E  is  the  electrostatic  held  between  the 
electrodes,  the  electrons  are  line-tied  and  can  only  Exit  drift  parallel  to  the  electrodes. 
In  principle  [18,  19],  this  prevents  the  electrons  from  screening  out  the  electric  held  in 
the  bulk  of  the  plasma  between  the  electrodes.  The  ions,  being  unmagnetized,  are  then 
accelerated  in  the  unscreened  electric  held.  This  simple  description  assumes  that  electrons 
are  collisionless.  If  there  is  a  finite  electron-ion  and/or  electron-neutral  collision  frequency, 
the  electrons  can  diffuse  across  held  lines.  Including  both  ion  and  electron  collisionality  in  a 
ID  fluid  theory  analysis,  Keidar  et  al.  argue  that  a  steady-state  spatial  distribution  results 
in  which  the  plasma  density  is  reduced  near  the  anode  side.  The  amount  of  density  reduction 
will  be  a  function  of  the  plasma  properties,  electrode  voltage,  and  magnetic  held  strength. 
The  aperture  of  the  antenna  used  for  communications  would  then  be  placed  directly  below 
the  reduced  density  region,  where  the  plasma  frequency  will  be  lower.  A  communication 
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FIG.  4:  ReComm  communication  through  plasma  layer  using  applied  electric  and  magnetic  fields. 
Adapted  from  Ref.  [17]. 

frequency  of  1  GHz  corresponds  to  a  free-space  wavelength  of  30  cm.  The  antenna  aperture 
will  be  on  the  order  of  the  wavelength,  w'hich  means,  as  a  rough  estimate,  we  can  assume 
that  the  electrode  separation  L  must  be  on  the  order  of  a  meter. 

In  this  report  we  describe  the  results  of  ID  simulations  of  cross-field  parallel-plate  simula¬ 
tions  in  the  plasma-optic  regime.  We  consider  only  the  issue  of  the  temporal  evolution  of  the 
plasma  density  in  the  region  between  the  electrodes  due  to  the  applied  voltage  and  magnetic 
field.  The  details  of  electromagnetic  wave  propagation  through  the  finite-thickness  plasma 
sheath  are  not  considered.  2D  and  3D  effects  are  neglected  for  the  present.  In  Sec.  Ill  A 
we  describe  results  from  a  baseline  ReComm  simulation  using  parameters  believed  to  be 
consistent  with  those  in  the  radio  blackout  regime.  The  simulation  results  suggest  that  the 
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Plate  separation 

L  (cm) 

1.0 

Applied  voltage 

E(V) 

100 

Applied  magnetic  field 

D{  G) 

2000 

Initial  plasma  density 

n  (cm-3) 

1011 

Initial  Electron  temperature 

Te  (eV) 

10.0 

Initial  Ion  temperature 

Ti  (eV) 

0.1 

Ion-neutral  collision  frequency 

Gn  (ns-1) 

0.066 

Electron-neutral  collision  frequency 

Ven  (ns-1) 

0.329 

TABLE  I:  Physical  parameters  for  the  baseline  ID  ReComm  simulation.  A  spatially  uniform  singly 
ionized  Ar  plasma  is  placed  between  plate  electrodes.  Results  from  the  simulation  are  shown  in 
Fig.  5. 

system  behaves  roughly  like  an  ambipolar  diffusion  process,  and  that  the  temporal  evolution 
of  the  plasma  occurs  on  corresponding  time  scales.  In  Sec.  Ill  B  we  describe  the  results  of 
a  series  of  parametric  studies  in  which  diode  length,  voltage,  and  plasma  collisionality  are 
varied  from  the  baseline  simulation  values. 


A.  ID  Simulations 

The  simulations  consist  of  a  ID  grid  with  conductors  at  the  boundaries.  A  constant 
voltage  difference  (after  a  short  temporal  ramp)  of  V  is  applied  between  the  cathode  (.r  =  0) 
and  anode  (x  =  L).  Table  I  shows  the  physical  parameters  used  in  what  we  will  refer  to  as 
the  “baseline”  ReComm  simulation.  In  this  section  we  will  discuss  results  from  this  single 
baseline  simulation.  Results  from  a  series  of  parametric  studies,  in  which  various  parameters 
are  varied  from  the  baseline  values,  will  be  discussed  in  the  following  section. 

An  initially-uniform  singly-ionized  Ar  plasma  is  placed  between  the  plates.  The  plasma 
density  of  1011  cm-3  corresponds  to  an  electron  plasma  frequency  of  u>pe  ~  20  GHz  which 
is  much  larger  than  the  L-band  frequencies  available  for  T&E  purposes.  We  assume  that 
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the  plasma  is  initially  at  rest  with  respect  to  the  electrodes.  That  is,  we  assume  that 
plasma  physics  phenomena  of  interest  occur  on  time  scales  short  compared  to  the  ratio  of 
the  plate  separation  and  the  hypersonic  How  velocity.  We  will  consider  the  validity  of  this 
assumption  later  in  the  report.  Ion-neutral  and  electron-neutral  collisions  are  included  in 
the  PIC  simulations  by  a  Monte-Carlo  collision  model  wdiich  has  been  described  in  detail 
elsewhere[20  -22],  For  simplicity  we  have  assumed  a  constant  ion-neutral  (i/in)  and  electron- 
neutral  (ven)  collision  frequency  (the  collision  frequency  should,  of  course,  in  general  be  a 
function  of  the  particle  energy).  Ionization  and  other  more  complicated  gas  chemistry  effects 
are  neglected  for  the  present.  For  this  reason,  it  is  not  very  significant  that  we  have  replaced 
the  N2  and  O2  molecular  ions  species  which  comprise  the  bulk  of  an  air  plasma  with  Ar  ions, 
which  are  only  slightly  heavier.  The  electron-collision  frequency  value  used  in  the  baseline 
simulation  corresponds  to  a  value  of  ven  ~  0.02icpe,  wdiich  is  in  accord  with  estimates  of  the 
value  in  the  hypersonic  plasma  sheath  [6,  12].  Since  wc  have  not  found  any  experimental 
data  or  estimates  for  the  ion-neutral  collision  frequency,  wc  have  chosen  to  set  vln  —  u)pi. 
where  tep,  is  the  ion  plasma  frequency.  Wc  note  that  a  collision  frequency  of  this  order  was 
found  to  damp  ion-acoustic  waves  driven  by  sheath  oscillations  near  the  cathode  and  anode 
surfaces  in  initial  simulations  with  negligible  electron  collisionality. 

For  these  initial  ID  simulations  we  wanted  to  attempt  to  resolve  electron  time  and  length 
scales  while  also  allowing  the  simulations  to  run  out  long  enough  to  observe  ion  dynamics. 
For  this  reason  we  began  with  a  relatively  small  plate  separation  value  of  L  =  1  cm.  If  we 
estimate  the  electric  field  in  the  gap  by  V/L  and  apply  a  2000  G  magnetic  field  as  shown 
in  Table  I,  the  E  x  B  drift  velocity  is  vp  ~  5  x  104  m/s,  which  is  about  10  times  larger 
than  the  ion  thermal  velocity  vti  and  10  times  smaller  than  the  electron  thermal  velocity  vte 
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(vti  <  Vd  <  vte)-  For  each  plasma  species  s  (=  e  or  i)  we  define  the  Larinor  radius  as 

max  [vD,  vts 
tls  = - 1 - 

where  ljcs  is  the  cyclotron  frequency  of  species  s.  We  hnd  that  for  the  parameters  in  the 
table  ~  0.04  cm  and  ru  ~  10  cm,  so  these  parameters  are  seen  to  be  in  the  plasma-optic 
regime  for  a  plate  separation  of  1  cm.  We  note  that  for  this  field  value  the  electrons  will 
still  make  many  cyclotron  orbits  in  between  collisions,  as  ven/u)ce  ~  0.01.  By  contrast  the 
ions  are  in  the  opposite  regime:  ~  100. 

By  choosing  a  plate  separation  of  1  cm,  we  can  choose  a  cell  size  of  Ax  =  0.04  cm. 
which  is  on  the  order  of  the  electron  cyclotron  radius,  without  requiring  a  prohibitively 
large  number  of  computational  grid  cells.  Several  hundred  plasma  particles  per  cell  were 
used  in  these  simulations  to  provide  a  good  statistical  representation  of  the  electron  and 
ion  energy  distribution  functions.  Macro  particles  which  strike  the  conducting  boundaries 
are  removed  from  the  simulation  without  replacement.  A  time  step  of  At  =  0.025  ns  gives 
uiceAt  ~  2u>peAt  ~  0.9,  ensuring  numerical  stability.  Simulations  with  this  time  step  value 
agree  well  with  simulations  at  smaller  time  step  values,  however,  a  significantly  larger  time 
step  would  require  the  use  of  implicit  algorithms  (Ref.  [23]  describes  a  new  algorithm  which 
allows  numerically  stable  simultaneous  under- resolution  of  uip  and  ujc).  This  would  probably 
be  required  for  larger  scale  simulations  (i.e.,  2D)  than  those  presented  here.  The  shortest 
ion  timescale  is  given  by  a ~  15  ns,  so  the  simulations  must  be  run  for  many  time  steps 
to  follow  the  ion  dynamics. 

Figure  5  shows  plots  of  ion  density  (top)  and  potential  (bottom)  across  the  anode-cathode 
gap  from  the  baseline  simulation  at  several  time  values.  At  t  =  50  ns  (a  few  ion  plasma 
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periods),  the  ion  density  is  uniform  except  at  the  electrode  surfaces  where  the  usual  Debye 
sheaths  have  been  established.  The  potential  profile  at  50  ns  shows  that  electric  field  across 
the  bulk  of  the  plasma  is  fairly  constant  and  ~  V j L.  Locally,  near  the  anode  and  cathode 
the  potential  decreases  as  the  surface  is  approached,  so  the  electrons  are  held  back  by  local 
sheath  electric  Helds.  In  general  the  sheath  electric  Helds  adjust  themselves  to  repel  the 
more  mobile  species.  In  this  way  quasi-neutrality  can  be  maintained  in  the  plasma.  The 
potential  drops  near  the  conducting  surfaces  are  of  the  order  of  the  electron  temperature. 

As  time  increases,  the  ion  density  profile  slowly  evolves  to  a  roughly  sinusoidal  shape, 
although  the  position  of  maximum  ion  density  shifts  towards  the  cathode  in  this  case.  So 
there  is  a  shifting  of  the  center-of-inass  of  the  plasma  towards  the  cathode  side.  However, 
significant  changes  in  the  ion  distribution  occur  only  on  time  scales  of  ~  10  /cs.  This  requires 
the  simulation  to  be  runs  for  ~  106  time  steps.  At  these  much  later  times  there  is  no  longer 
a  (roughly)  constant  electric  Held  in  the  bulk  plasma.  There  is  a  fairly  small  electric  Held 
on  the  cathode  side  of  the  simulation  with  a  larger  electric  Held  (>  V/L)  near  the  anode. 
We  note  that,  other  than  the  slight  skewing  of  the  density  distribution  toward  the  cathode, 
the  general  temporal  behavior  resembles  that  of  a  simple  ID  ambipolar  diffusion  problem 
between  parallel  plates  without  an  applied  voltage. 

This  situation  is  discussed  at  some  length  in  Ref.  [24],  In  Appendix  B  we  extend  the 
analysis  of  Ref.  [24]  to  include  the  effects  of  an  applied  voltage.  In  the  theory  outlined  in 
Appendix  B  the  total  plasma  number  density  is  assumed  to  decrease  exponentially  with  a 
time  constant  of  r  for  the  lowest  order  diffusion  mode.  We  show  in  the  following  section, 
that  this  is  approximately  the  case  for  the  simulations  presented  in  this  report.  For  given 
geometry  and  plasma  properties,  the  time  constant  r  can  be  calculated  as  a  function  of 
the  three  voltage  scales,  the  applied  voltage  V,  the  ambipolar  potential  Va  [Eqs.  (B5)  and 
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FIG.  5:  Snapshots  of  ion  density  (top)  and  potential  (bottom)  as  a  function  of  x  at  several  times 
for  the  baseline  ID  ReComm  simulation  with  parameters  given  in  Table  I. 

(B7)],  and  Vm,  the  potential  difference  between  the  center  of  the  diode  and  the  cathode 
[Eq.  (B13)].  The  applied  voltage  V  is,  of  course,  a  known  input  parameter.  The  ambipolar 
potential  Va  is  uniquely  determined  by  the  bulk  plasma  properties  and  is  assumed  to  be 
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constant  in  the  simple  diffusion  theory,  and  so  can  in  principal  be  calculated  from  the  initial 
conditions  of  the  plasma.  But  the  mid-potential  Vm  cannot  be  determined  from  a  quasi¬ 
neutral  diffusion  theory  and  depends  on  the  details  of  sheath  physics  at  electrode  surfaces. 
In  the  usual  treatment  of  ambipolar  diffusion  it  is  approximated  that  |  VX/ /V^|  1  (this  is 

equivalent  to  the  statement  that  the  plasma  density  at  the  electrodes  is  vanishingly  small 
at  the  electrodes),  the  timescale  r  is  given  approximately  by 


r  = 


(52) 


where  Da  is  the  ambipolar  diffusion  coefficient  defined  by  Eqs.  (B6)  and  (B7)  [see  Appendix 
A].  If  we  calculate  r  using  the  parameters  in  Table  I,  we  obtain  r  ~  50/is.  This  simple 
estimate  yields  a  timescale  ~  10  /is,  wrhieh  is  in  rough  accord  with  the  simulation  results. 
Again  using  the  parameters  of  the  table,  we  can  also  estimate  Va  ~  — 5.G  eV.  The  analysis 
in  Appendix  B  demonstrates  that  when  an  applied  voltage  is  present  and  Va  <  0  (meaning 
that  the  electrons  are  more  mobile),  the  center  of  mass  of  the  density  distribution  should 
shift  towards  the  cathode.  This  is  also  observed  in  the  baseline  simulation  (Fig.  5). 

We  note  that  the  diffusion  theory  outlined  in  Appendix  B  makes  the  assumption  that 
both  Va  and  Da  are  independent  of  both  position  and  time.  This  implies  that  electron  and 
ion  temperatures  and  collision  frequencies  are  constant  as  well.  We  consider  nowr  the  time 
variation  of  the  temperature.  Figure  6  shows  the  electron  and  ion  temperature  (in  units  of 
volts)  as  well  as  Vm  and  Va  (calculated  by  using  instantaneous  values  of  temperature)  as  a 
function  of  time.  We  note  that  the  temperature  histories  shown  here  represent  a  measure  of 
velocity  spread  over  the  entire  simulation  space.  Any  spatial  variations  in  temperature  are 
averaged  out  in  this  process. 
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FIG.  6:  Spatially  averaged  electron  and  ion  temperatures  and  voltage  scales  as  a  function  of  t  for 
the  baseline  ID  ReComm  simulation  with  parameters  given  in  Table  I.  The  ambipolar  voltage  scale 
Va  is  defined  by  Eq.  (B5),  and  Vm  is  the  potential  difference  between  the  center  and  the  cathode 
[Eq.  (B13)]. 

From  Eq.  (B8)  we  see  that  Va  is  a  linear  combination  of  Te  and  Tt  (for  fixed  magnetic  field 
and  collision  frequencies).  Note  that  Va  can  be  positive  or  negative  depending  on  the  relative 
mobility  of  ions  and  electrons.  The  simple  diffusion  theory  outlined  in  the  appendix,  assumes 
that  both  Va  and  the  diffusion  coefficient  D„  are  constant  in  space  and  time.  This  implies 
that  plasma  temperatures  and  collision  frequencies  are  constant.  While  collision  frequencies 
are  fixed  by  the  code,  we  see  from  Fig.  6  that  the  average  simulation  temperatures  exhibit 
some  heating  until  about  10  /is,  after  wffiich  they  level  off.  We  restrict  our  consideration  to 
times  greater  than  10  /is  and  neglect  the  time  dependence. 

When  using  the  initial  temperatures,  we  find  that  Va  ~  —5.6  eV.  If  we  use  instead  the 
asymptotic  temperatures,  Te  ~  15  eV  and  T,  ~  6  eV,  the  ambipolar  potential  changes  only 
slightly  Va  ~  —6  eV.  From  the  simulation  we  have  an  asymptotic  value  of  Vm  ~  Va.  As 
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shown  in  Appendix  B.  for  |V/\4|  3>  1.  a  physical  solution  to  the  diffusion  equation  requires 
that  VM/Va  <  ln{ 2).  From  Fig.  6,  this  clearly  not  the  case. 

This  apparent  contradiction  can  be  explained  as  follows.  We  again  emphasize  that  in 
Fig.  6  Va  is  calculated  by  using  spatially  averaged  values  of  the  species  temperatures.  Since 
the  electric  fields  are  not  uniform  across  the  gap,  there  is  actually  an  x  variation  in  tem¬ 
perature  due  to  non-uniform  Joule  heating.  This  can  be  seen  in  Fig.  7,  where  the  bottom 
plot  shows  the  electron  and  ion  temperatures  (in  units  of  volts)  as  a  function  of  x  at  t  =  20 
//s.  We  see  that  electrons  are  heated  fairly  uniformly  across  the  gap,  but  ions  are  heated 
preferentially  on  the  anode  side,  where  the  electric  held  is  largest.  We  have  also  plotted 
Va,  calculated  by  Eq.  (BG)  using  the  local  values  of  temperature,  as  a  function  of  x.  The 
spatial  variation  in  Va  (and  Da)  implies  that  the  system  no  longer  satisfies  a  simple  diffusion 
equation. 

The  upper  plot  of  Fig.  7  shows  both  the  electron  and  ion  density  distributions  at  t  =  10 
/is.  The  data,  which  come  from  an  instantaneous  snapshot,  are  somew'hat  noisy.  But  it  is 
clear  that  there  is  a  significant  charge  imbalance  occurring  in  a  region  which  extends  from 
the  anode  about  2  mm  into  the  bulk  plasma.  This  is  therefore  a  large  region  (much  larger 
than  a  Debye  length)  where  the  assumption  of  charge  neutrality  does  not  hold.  In  this  same 
region  Va  is  positive,  which  implies  that,  at  least  locally,  the  ions  are  more  mobile  than  the 
electrons.  From  ambipolar  theory  this  would  predict  an  ambipolar  field  in  this  region  which 
would  hold  back  the  electrons  from  the  anode,  which  is  clearly  not  the  case  (as  seen  from 
the  potential  plot).  But  the  large  electric  held  is  not  ambipolar,  since  there  is  a  significant 
charge  imbalance.  Alternately,  in  this  region  electron  and  ion  fluxes  are  not  nearly  equal. 
This  flux  imbalance  is  required  to  draw  current  from  the  cathode  to  the  anode. 

As  discussed  in  Appendix  B,  the  simple  diffusion  theory  also  predicts  a  constant  ion 
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FIG.  7:  Results  from  the  baseline  simulation  at  1  —  20  //s.  (Top)  Ion  and  Electron  density  profiles. 
(Bottom)  Ion  and  Electron  temperature,  Va,  and  the  electrostatic  potential  4>  (all  in  units  of  volts) 
as  a  function  of  x. 

center-of-mass  once  the  higher-order  modes  have  decayed,  although  the  density  profiles  in 
Fig.  5  exhibit  a  cont  inuing  drift  of  the  center-of-mass  at  late  times.  This  can  also  be  explained 
by  the  fact  that  Va  and  Da  are  not  spatially  uniform.  In  this  case  it  is  no  longer  possible 
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to  exactly  separate  the  variables  x  and  t,  of  the  density  function,  as  in  done  Eq.  (B9).  The 
combination  of  this  effect  and  the  breakdown  of  quasi-neutrality  leads  to  a  time  dependence 
in  the  position  of  the  plasma  center-of-mass.  This  is  seen  more  clearly  in  the  simulation 
results  presented  in  Sec.  Ill  B. 

We  have  demonstrated  here  that  qualitatively  our  baseline  ReComm  simulation  behaves 
as  an  ambipolar  diffusion  process,  in  which  the  density  profile  evolves  on  time  scales  con¬ 
sistent  with  diffusion  theory.  The  density  profile  evolves  to  a  roughly  sinusoidal  shape,  in 
which  the  center-of-mass  is  shifted  due  to  the  applied  voltage.  The  existence  and  sign  of 
this  shift  can  also  be  predicted  by  the  simple  diffusion  theory.  As  we  have  demonstrated, 
however,  there  is  significant  Ohmic  heating  and  large-scale  charge  separation  present  in  the 
simulations,  which  preclude  using  the  theory  quantitatively.  In  the  following  sections  we 
present  the  results  of  a  series  of  parametric  studies  in  which  the  input  parameters  of  the 
ReComm  simulation  are  varied  systematically. 

B.  Parametric  Studies 

In  this  section  we  describe  the  results  of  parametric  studies  in  w'hich  certain  physical 
parameters  are  varied  from  the  baseline  values  of  Table  I.  We  begin  by  describing  the 
results  from  a  series  of  simulations  in  which  the  electrode  separation  L  is  varied.  All  other 
physical  parameters  remain  unaltered  from  the  values  in  the  table.  We  pointed  out  in  the 
previous  section  that  on  long  time  scales  the  time  evolution  of  the  plasma  density  distribution 
resembles  an  ambipolar  diffusion  process.  The  timescale  for  a  standard  ambipolar  diffusion 
process  (he.,  with  no  applied  field)  should  scale  as  L2 / Da.  Figure  8  shows  results  for 
simulations  with  L  =  0.5, 1. 1.5, 2, 3,  and  4  cm.  As  mentioned  above,  the  ion  Larmor  radius 
for  the  present  set  of  input  plasma  parameters  and  magnetic  field  is  ~  10  cm.  Simulations 
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with  larger  plate  separations  would  require  a  decrease  of  the  magnetic  held  strength  to 
remain  in  the  plasma-optic  regime.  Since  the  cell  size  Ax  has  been  maintained  for  all  these 
simulations,  runs  at  larger  values  of  L  are  also  more  time  consuming.  The  top  plot  shows  the 
total  simulation  ion  charge  Q(t)  [normalized  by  the  initial  charge  value  Q(0)]  as  a  function  of 
time.  For  all  values  of  L,  after  a  transient  time,  the  total  ion  charge  decreases  exponentially. 
This  is  consistent  with  diffusion  of  the  lowest  order  mode. 

In  the  diffusion  model,  when  higher  order  diffusion  modes  have  decayed  away,  Q(t)  is 
proportional  to  the  integral  over  x  of  n(x,t),  as  given  by  Eq.  (B9)  of  Appendix  B.  So  after 
the  transient  time,  the  total  ion  charge  can  be  expressed  as 

Q(t)  oc  e~t/Tf.  (53) 

For  each  simulation,  the  timescale  tj  is  found  numerically  by  fitting  the  results  for  Q(t). 
An  alternate  timescale  can  also  be  calculated  from  the  first  e-folding  time  of  the  total  ion 
charge: 

Q(re)/Q(0)  =  e~\  (54) 

where  Q( 0)  is  the  initial  ion  charge.  Both  Tj  and  re  are  plotted  as  a  function  of  L  (the 
anode-cathode  gap  spacing)  in  the  bottom  plot.  For  the  smaller  values  of  L  both  r  values 
scale  approximately  as  L 2 ,  but  for  L  >  2  cm,  the  quadratic  scaling  falls  off  to  something 
closer  to  a  linear  dependence.  This  can  be  explained  by  referring  to  ambipolar  diffusion 
theory.  As  shown  in  Appendix  B.  the  dimensionless  parameter  L2  /  Dar  is  a  function  of  both 
V )Va  and  V\ t/Va.  Da  and  Va  are  functions  only  of  the  bulk  plasma  parameters,  however, 
the  mid  potential  Vm,  which  cannot  be  determined  bv  simple  quasi-neutral  fluid  theory,  is 
in  general  a  function  of  L.  This  L  dependence,  wdiich  can  be  seen  in  the  bottom  plot  of 
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Fig.  9,  breaks  the  simple  scaling  relation  r  oc  L2. 

Figure  9  shows  the  ion  density  and  potential  for  each  simulation  as  a  function  of  xj L  at 
the  time  t  =  re.  The  density  profiles  have  a  nearly  self-similar  shape,  except  for  a  slight 
shifting  of  the  ion  center-of-mass  back  toward  the  center  of  the  diode  with  increasing  L. 
The  potential  profiles  have  qualitatively  the  same  shape,  but  the  mid  potential  Vm  increases 
somewhat  with  increasing  L.  We  also  note  that  for  each  run  the  mid  potential  Vm  is  roughly 
constant  in  time  after  a  transient  time  of  order  re.  The  same  is  also  roughly  true  of  the 
average  electron  and  ion  temperatures. 

Figure  10  shows  ion  (x) / L  and  xrms/L  as  a  function  of  f/re  for  all  the  simulations  with  L 
varied.  For  all  runs  the  ion  center-of-mass  is  skewed  towards  the  cathode.  As  noted  in  the 
previous  section,  for  a  given  L ,  the  center-of-mass  drifts  monotonically  towards  the  cathode 
in  time.  At  a  hxed  value  of  normalized  time  t/re  the  ion  center-of-mass  shifts  back  towards 
the  center  of  the  two  plates  as  L  is  increased.  So  the  skewing  of  the  density  distribution 
diminishes  as  L  is  increased  (for  a  fixed  magnetic  held  strength).  The  plot  xrms/L  shows 
that  the  “spread”  in  ion  density  is  not  a  strong  function  of  L.  This  can  also  be  clearly  seen 
in  the  top  plot  of  Fig.  9.  Note  also  that  unlike  the  center-of-mass,  the  density  spread  does 
seem  to  asymptote  for  t  >  re. 

We  now  consider  a  series  of  runs  for  which  only  the  applied  voltage  V  is  varied.  All  other 
simulation  parameters  are  the  same  as  in  Table  I,  including  the  electrode  plate  separation 
which  is  hxed  at  1  cm.  The  top  plot  of  Fig.  11  shows  results  for  normalized  total  ion  charge 
as  a  function  of  time  for  runs  with  V  =  0,  50, 100,  200  V.  Note  that  for  the  V  =  0  case, 
which  corresponds  directly  to  the  ambipolar  diffusion  problem  considered  in  Ref.  [24],  the 
time  dependence  of  the  total  charge  is  not  purely  exponential  at  late  times.  This  is  due  to 
temperature  gradients  and  temporal  variations  in  both  V\ /  and  the  temperatures. 
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FIG.  8:  Results  from  a  series  of  ID  ReComm  simulations.  The  simulation  parameters  are  the  same 
as  those  given  in  Table  I,  except  that  the  plate  separation  L  is  varied.  (Top)  Normalized  total  ion 
charge  as  a  function  of  time  for  various  values  of  L.  (Bottom)  Time  scales  tj  and  re  as  a  function 
of  L. 
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FIG.  9:  Results  from  a  series  of  ID  ReComm  simulations.  The  simulation  parameters  are  tile 
same  as  those  given  in  Table  I,  except  that  the  plate  separation  L  is  varied.  (Top)  Normalized  ion 
density  at  t  =  re  as  a  function  of  x/L.  (Bottom)  Potential  at  t  =  re  as  a  function  of  x/L. 

For  each  run,  the  values  of  tj  and  re  are  calculated  from  the  Q(t)  data.  Of  course  the 
value  of  tj  for  V  =  0  is  a  rough  estimate  by  necessity.  The  inverse  values  of  re  and  17  are 
plotted  as  a  function  of  V  in  the  bottom  plot  of  Fig.  11.  For  V  =  50, 100,  and  200  V,  Va  <  0 
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FIG.  10:  Results  from  a  series  of  ID  ReComni  simulations.  The  simulation  parameters  are  the 
same  as  those  given  in  Table  I,  except  that  the  plate  separation  L  is  varied.  (Top)  {x)/L  as  a 
function  of  t/re.  (Bottom)  xTms/L  as  a  function  of  t / re . 

and  V  »  1141.  In  this  regime,  the  data  suggests  that  r  is  roughly  proportional  \/V .  That  is, 
the  plasma  density  can  be  moved  around  more  quickly  as  the  electrode  voltage  is  increased. 
Figure  12  shows  the  ion  density  and  potential  for  each  simulation  as  a  function  of  x  at 
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FIG.  11:  Results  from  a  series  of  ID  ReComm  simulations.  The  simulation  parameters  are  the 
same  as  those  given  in  Table  I,  except  that  the  plate  voltage  V  is  varied.  (Top)  Normalized  total 
ion  charge  as  a  function  of  time  for  various  values  of  V.  (Bottom)  Time  scales  ry  and  re  as  a 
function  of  V. 
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the  time  t  =  re.  The  density  profiles  again  have  a  nearly  self-similar  shape,  except  for  the 
shifting  of  the  ion  center-of-mass  toward  the  plate  for  the  runs  with  an  applied  voltage. 
With  no  applied  voltage,  the  distribution  is  symmetric  about  the  midpoint  between  the 
plates,  as  is  expected  by  symmetry.  If  we  neglect  transients,  the  position  dependence  of  the 
temperatures,  as  well  as  charge-separation  effects,  we  can  roughly  estimate  that  Va  ~  —2.5  V 
and  Vjvf  ~  10  V  for  the  run  with  V  =  0.  For  these  parameters  we  can  estimate  [by  solving 
Eq.  (B18)]  that  r  ~  300  /is.  This  is  in  rough  agreement  with  the  measured  value  of  Tj.  For 
this  case  (V  =  0),  we  note  that  the  mid-potential  Vm  is  of  the  order  of  Te/e  which  is  a  few 
times  —  Va.  The  diffusion  theory  of  Ref.  [24]  predicts  an  infinite  value  of  Vm  for  a  boundary 
condition  with  vanishing  plasma  density  at  the  electrode  surfaces.  This  does  not  actually 
occur  as  the  assumption  of  quasi-neutrality  breaks  down  near  the  electrodes  and  limits  the 
magnitude  of  Vm-  Figure  13  shows  the  plasma  center-of-mass  (.r)  and  spread  xrms  for  each 
run  as  a  function  of  //re.  The  temporal  drift  of  the  center-of-mass  is  again  evident,  except, 
of  course,  for  V  =  0  case. 

In  the  next  series  of  runs,  the  electron-neutral  collision  frequency  ven  is  varied  from  its 
baseline  value  of  5 ujpi.  We  have  performed  simulations  for  i'en/uPi  =  1,2, 3.5, 5,  and  10.  The 
calculated  time  constants  re  and  tj  are  plotted  as  a  function  of  uen/uipi.  We  note  that  with 
increasing  electron  collisionality,  the  time  constants  begin  to  level  off.  As  ven  is  increased,  the 
electrons  become  more  and  more  mobile  [for  a  magnetized  plasma  species  w'here  locs  /pTI 
the  diffusion  coefficient  and  mobility  both  scale  linearly  with  i/sn;  see  Eq.  (B7)].  In  the  limit 
that  the  electrons  become  much  more  mobile  than  the  ions,  ambipolar  diffusion  becomes 
limited  by  ion  mobility.  For  this  reason  r  oc  1  / Da  asymptotes  as  seen  in  Fig.  14. 

Figure  15  shows  the  ion  density  and  potential  as  a  function  of  x  at  t  =  re  for  each  of  the 
simulations.  We  note  first  that  for  the  smaller  values  of  ven/u>pi ,  for  which  the  electrons  are 
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FIG.  12:  Results  from  a  series  of  ID  ReComm  simulations.  The  simulation  parameters  are  the 
same  as  those  given  in  Table  I.  except  that  the  plate  voltage  V  is  varied.  (Top)  Normalized  ion 
density  at  t  =  Te  as  a  function  of  x.  (Bottom)  Potential  at  t  =  Te  as  a  function  of  x. 

relatively  less  mobile,  there  large  are  regions  near  the  cathode  completely  devoid  of  ions. 
Since  the  electron  density  is  not  negligible  in  this  region,  there  is  again  a  large  electrode 
sheath  where  charge  neutrality  breaks  down.  For  these  runs  the  potential  plots  show  that. 
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FIG.  13:  Results  from  a  series  of  ID  ReComm  simulations.  The  simulation  parameters  are  the 
same  those  as  given  in  Table  I,  except  that  the  plate  voltage  V  is  varied.  (Top)  (x)  as  a  function 
of  t/re.  (Bottom)  xrms  as  a  function  of  £/re. 

most  of  the  potential  is  screened  out  in  the  region  near  the  cathode,  with  only  a  small  electric 
in  the  remainder  of  the  gap.  In  the  limit  that  the  ven  — >  0,  the  electrons  are  completely 
immobile,  as  they  are  line-tied  and  can  only  drift  perpendicular  to  the  electrodes.  In  this 
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FIG.  14:  Timescales  ry  and  re  are  plotted  as  a  function  of  venju)v%.  The  simulation  parameters 
are  the  same  as  those  given  in  Table  I,  except  that  the  electron-neutral  collision  frequency  vcn  is 
varied. 

case,  the  potential  is  completely  screened  out  in  a  region  near  the  cathode  from  which  ions 
are  excluded.  In  the  remainder  of  the  gap  there  is  no  electric  field  and  the  ion  density 
distribution  remains  unaltered.  As  the  electron  collision  frequency  is  increased  relative  to 
the  ions,  the  large  cathode  sheath  eventually  recedes. 

Figure  16  shows  the  ion  density  center-of-mass  and  spread  as  a  function  of  t/re.  Note 
that  the  center-of-mass  shifts  towards  the  anode  for  low  electron  collisionality,  but  shifts 
back  towards  the  cathode  when  electron  collisionality  is  increased.  This  agrees  with  the 
general  rule,  established  Appendix  B,  that  the  direction  of  the  shift  of  the  center-of-mass 
is  determined  by  the  relative  mobility  of  the  two  species.  When  the  electrons  are  more 
mobile,  the  density  shifts  towards  the  cathode,  and  when  the  ions  are  more  mobile  the  shift 
is  towards  the  anode. 
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FIG.  15:  Results  from  a  series  of  ID  ReCoinm  simulations.  The  simulation  parameters  are  the 
same  as  those  given  in  Table  I,  except  that  the  electron-neutral  collision  frequency  vtn  is  varied. 
(Top)  Normalized  ion  density  at  t  =  re  as  a  function  of  x.  (Bottom)  Potential  at  t  =  re  as  a 
function  of  x. 
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where 
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J1  —  u;2 


(58) 


and 


u): 


rj=  l-~f. 

u r 


(59) 


The  Green’s  function  for  Eq.  (57)  is 


G(f)  = 


-2  [\er]\  (7 V  -  x2  -  y2)]  1/2  ,  |7z|  >  yG:2  +  y2, 
2i  [|f77|  (7  V  -  x2  -  y2)]-1/2  ,  I72I  <  v^+y2. 


(60) 


where  72  =  —e/77. 

For  the  charge  density  in  an  antenna  defined  as 


70cos[k(L  -  |x|)] 

p{x,y,z)  =  — - ,  - sgn  (x)6(y)6(z), 

c  sm  (kL) 


(61) 


we  can  obtain  the  following  expression  for  the  electrostatic  potential  of  oscillations  excited 
by  the  dipole  antenna  inside  a  plasma  sheath: 


<t>(r)  = 


2_  fL _ p(Q<% _ 

fv\  Jl  \J 7222  -  {x  -  O2  -  y2 


(62) 


The  electrostatic  potential  defined  by  Eq.  (62)  contains  resonances  surfaces  where  it  has 
singularities.  To  find  these  surfaces  and  the  type  of  electrostatic  potential  singularity  on 
these  surfaces,  we  define  the  function 


=  7 V  -  {x  -02  -  y2- 


(63) 
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On  the  resonance  surface,  the  electrostatic  potential  $  a  logarithmic  singularity.  This  sin¬ 
gularity  can  be  identified  from  the  conditions 


M£,x,y,z)  =  0, 


(64) 


dfi 

dZ 


-2(*-O  =  0. 


(65) 


These  conditions  (64)-(65)  define  parts  of  the  plane 


72  =  ±y, 


(66) 


with  —  L  <  x  <  L.  Now  f\  can  be  put  into  the  form 


/i(£>#>  y,  z)  = 


4  f\ (^si  xs,  ys,  2S)  ( e  .2 


d? 


(£.  -  xisY+Sf, 


(67) 


where 


X,  Bh 
Sf=  to 


(x-xs)+ 


df  1 


X~X8 


dy 


{y-ysY 


dfi 


y=y* 


dz 


{z-zs)  =  2'yzs{z—za)—2ys(y—ys).  (68) 


z—z» 


With  Eq.  (67),  we  can  now  write  the  following  expression  for  the  potential  of  the  quasi¬ 
electrostatic  whistler  wave  in  the  vicinity  of  the  resonance  surface: 


<*>,  =  - 


2pi 

\/N 


/ 


2~gjjpr  “  x?o)2  +  Sf 


1-1/2 


d£. 


(69) 
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Carrying  out  the  integration  in  Eq.  (69)  with  8f  ^  0,  we  find 


=  - 


-)2=l  n(Uf 

N  V2 


+  constant. 


(70) 


From  Eq.  (70)  it  is  clearly  seen  that  when  we  are  approaching  the  resonance  surface  and 
8f  Q,  the  potential  of  the  quasi-electrostatic  whistler  wave  experiences  a  logarithmic 
singularity.  We  can  use  Eq.  (70)  to  find  the  electric  held  excited  by  an  antenna  in  the 
vicinity  of  the  resonance  surface.  Introducing  the  coordinates  along  and  perpendicular  to 
the  resonance  surface, 


v'l  +  72 


and 


7  z-y 

v/TTfr 


(71) 


we  can  find  that  the  strongest  singularity  is  experienced  by  the  electric  field  component 
directed  the  normal  to  the  resonance  surface.  The  expression  for  the  electric  field  has  the 
following  form, 


k  cos[/c(L  -  |x|)]  , 

Ea  ~  21 A - 7= — —====  sgn(x), 


(72) 


W\/M  ^x/1  +72 

where  5  is  the  distance  from  the  resonance  surface.  From  the  expression  for  the  excited 
field,  Eq.  (72)  indicates  that  the  field  experiences  the  1/5  singularity  in  the  vicinity  of 
the  resonance  surface.  One  way  to  eliminate  this  singularity  is  to  include  the  finite  plasma 
temperature  into  the  components  of  the  dielectric  tensor.  Another  possibility  is  to  introduce 
the  finite  width  of  the  dipole  antenna.  Equation  (72)  is  valid  at  distances  from  the  resonance 
surface  that  exceed  the  electron  cyclotron  radius  or  finite  antenna  width  whatever  is  larger. 
This  analysis  allows  us  to  make  an  estimate  of  the  wave  field  excited  by  the  dipole  antenna 
near  the  resonance  surface. 


We  can  estimate  the  relative  importance  of  different  nonlinear  processes  in  vicinity  of 
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the  plasma  sheath.  Using  Eq.  (72)  as  an  estimate  of  the  electric  field  excited  by  the  dipole 
antenna  inside  the  plasma  sheath,  we  estimate  the  change  in  the  plasma  temperature  due  to 
the  presence  of  the  dipole  antenna  and,  eventually,  the  level  of  the  random  noise  produced, 
which  can  influence  the  quality  of  the  communication  channel. 

The  oscillation  frequency  in  the  antenna  is  comparable  to  the  lower  hybrid  resonance 
(LHR)  frequency,  Eq.  (55),  with  characteristic  wavelength  obtained  from  ( kc/ujpe )2  »  1. 
For  the  electron  density  inside  the  sheath,  ne  ~  1011  cm-3,  this  corresponds  to  the  excited 
LOR  oscillation  wavelengths  of  the  order  of  several  centimeters  (1  3  cm): 


A  <  A0  = 


2tcc 

Ldpe 


~  10  cm. 


(73) 


For  higher  plasma  densities,  these  exited  wavelengths  are  smaller.  The  excited  oscillations 
are  known  to  occur  inside  the  resonance  cones  oriented  along  the  external  magnetic  field 
with  an  apex  angle  a  ~  c o/uice. 

The  important  parameter  which  characterizes  the  nonlinear  processes  in  the  vicinity  of 
the  antenna  is  the  ratio  of  the  electric  field  energy  to  the  plasma  pressure.  For  the  values  of 
this  parameter  inside  the  range: 


777 


< 


E 2 


M  47T77oT 


<  1. 


(74) 


The  principal  nonlinear  effect  is  wave  energy  spectrum  migration  towards  smaller  wave¬ 
lengths  and  the  collapse  of  the  lower  hybrid  waves.  This  eventually  results  in  resonance 
absorption  of  the  wave  energy  by  electrons  and  ions.  When  the  excited  filed  is  stronger  and 
satisfies  the  condition, 


E2 

AirnoT 


(75) 
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strong  nonlinear  effects  will  be  under  way.  The  high  radiation  pressure  will  push  the  plasma 
out  from  the  resonance  cone  region  and  the  intersection  of  the  particle  trajectories  due  to  the 
strong  parametric  instability  will  heat  the  plasma  as  a  w'hole.  Taking  the  current  in  antenna 
I  a  ~  3  to  4  A,  and  using  the  plasma  parameters  and  magnetic  held  strength  presented 
above,  we  can  use  Eqs.  (72)  ami  (75)  to  show  that  the  LOR  oscillations  in  the  vicinity  of 
the  antenna  are  unstable. 

B.  Plasma  heating  by  quasi-electrostatic  whistler  waves  in  the  vicinity  of  a  VLF 
antenna  inside  a  plasma  sheath 

Assuming  that  the  wavelength  of  the  excited  waves  is  on  the  order  of  the  thickness  of  a 
plasma  sheath,  we  examine  the  nonlinear  interaction  of  the  pump  Held, 

E0(t)  =  E0s\n(uJot) ,  (76) 

with  the  plasma  sheath.  We  are  interested  in  the  case  when  the  waves  excited  by  the  antenna 
Held  is  strong  enough  to  satisfy  condition  (75),  indicating  the  presence  of  strongly  nonlinear 
coupling  effects.  The  high  radiation  pressure  will  push  out  the  plasma  from  the  resonance 
cone  region  and  the  intersection  of  particle  trajectories  due  to  a  strong  parametric  instability 
will  heat  the  plasma  as  a  whole.  The  instability  of  the  plasma  in  the  external  electric  and 
magnetic  Helds  due  to  excitation  of  electrostatic  oscillations  near  the  lower  hybrid  resonance 
was  investigated  in  a  series  of  papers.  The  dispersion  equation  for  parametric  instability 
for  the  case  when  the  oscillation  speed  of  electrons  in  the  external  electric  Held  is  small  in 
comparison  with  the  thermal  velocity,  including  the  influence  of  Coulomb  collisions,  was 
derived  by  Andreev  and  Kirii  [25].  Nonresonant  decay  neglecting  the  thermal  motion  of 
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particles,  which  is  possible  in  a  strong  electric  field  of  the  pmnp  wave  when  the  energy  of 
oscillations  of  electrons  is  sufficiently  larger  than  their  thermal  energy,  is  examined  by  Porko- 
lab  [26].  The  experimental  observation  of  parametric  instability  in  the  range  of  frequencies 
we  are  interested  in  and  with  the  frequency  of  the  pump  wave  uJq  <  3 Um  was  presented 
by  Bonizzoni  et  al.  [27],  where  it  was  noticed  that  the  excitation  of  a  broad  spectrum  of 
non-resonant  quasi-ion  modes  takes  place.  This  is  also  the  situation  in  our  case,  but  now 
with  a  very  strong  pump  wave,  w’hen  the  condition  (75)  is  satisfied  and  forced  oscillations 
with  the  growth  rate  7  (for  u;ri  <SC  7  <  oJm)  are  excited.  Contrary  to  the  case  considered 
by  Bonizzoni  et  al..  we  ignore  the  dispersive  effects  due  to  the  thermal  motion  of  plasma 
particles. 

C.  Dispersion  Equation  for  Parametric  Instability 

We  present  a  brief  derivation  of  a  dispersion  equation  which  describes  the  parametric 
excitation  of  lower  oblique  resonance  oscillations  in  the  presence  of  strong  pump  wave  with 
the  frequency  <-c0  >  oJlh-  We  consider  the  pump  field  homogeneous  in  space,  with  temporal 
profile  given  by  Eq.  (76).  As  mentioned  above,  the  wavelengths  under  consideration  are 
short  in  comparison  with  the  characteristic  scale  length  of  the  inhomogeneity  of  the  field  of 
the  pump  wave.  Taking  into  account  that  the  wavelength  of  the  quasi-electrostatic  whistler 
waves  excited  by  the  VLF  antenna  are  A  ~  27T c/Upe,  giving  k  cjpe/c.  We  analyze  nonlinear 
excitation  of  quasi-electrostatic  oscillations  with  the  frequencies  above  the  lower  hybrid 
frequency  due  to  a  presence  of  the  pump  field  produced  by  a  VLF  antenna.  Below',  ne  and 
n,  are  the  electron  and  ion  density  perturbations  in  the  excited  oscillations  and  Ve  and  Vj 
are  the  velocity  perturbations.  Furthermore,  Sne  and  5iii  as  well  as  S Ve  and  SVi  are  the 
perturbations  parametrically  connected  to  the  low  frequency  plasma  dynamics. 
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For  the  electron  speed  at  high  frequency  one  can  write: 


6  6/9  -*  -* 

VPT  —  Eu  —  (  ~  +  ikL  Oe  ) 


UlUJr, 


y  mu2  \dt 


(77) 


and 


Kv  =  - 


T  IVjJr 


■Ex--^T(~  +  ikU0e)Ey. 

rnuJt  \  dt  1 


(78) 


In  Eqs.  (77)  and  (78),  E  =  E(t)  exp(ikr)  is  the  parametrically  excited  field  of  oscillations. 
The  electron  speeds  due  to  the  pump  wave  are  defined  as 


U()e  —  (E()y€x  Eox&y)  sill  ((x?o^)  • 

77  UJpp 


(79) 


Now  we  can  write  expressions  for  the  perturbations  of  the  electron  and  ion  densities  due 
to  the  high  frequency  oscillations,  using  the  continuity  and  Poisson  equations, 


d_ 

dt 


u ’le  k2z 


1 


—  +  ikUoe  )  ne  +  7^nP  =  —  V  ,  o  (  +  ikUoe  )  n<+  (80) 


1  +  Lj2/u2e  k2  1  +  u2  /u)2e  \dt. 


d 


B0  1  +<«£/ u2e  \dt 


d  -Tft 
—  +  ikUoe 


Eoi  x  V5n 


wfe  m  /  d6n ) 

+  —  [ E0z  cos  (oj0t)}  —  >  , 


J  z  (Jo 


n 

—  +  ikUoi  )  nt=Jk  (ne  -  nt) . 


(81) 


In  Eq.  (81),  the  expression 


ft  eE0  , 

Uoi  =  - - COS  ((Jot) 

M  (Jo 


is  the  oscillatory  speed  of  the  ions  in  the  pump  wave  field.  We  now  introduce  the  function 


T]e(t)  =  ne  exp  (—ikSe'j  , 


(82) 
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which  corresponds  to  the  electron  charge  density  in  the  coordinate  frame  oscillating  together 
with  the  electrons  under  the  action  of  the  pump  held.  In  Eq.  (82),  for  ae  we  have 


ae(t)  = - {Eoyex  -  E0xey)  cos  (u0t)  + 

rriidceLJo 


- 2  E0z  sin  {u)0t)ez. 

juuJq 


(83) 


For  the  ions,  we  consider 


d  rrT 

-gt  >  kUoi. 


After  averaging  over  the  fast  time  ~w0',  we  obtain 


1  ,  ^ pe  \  9  He  ,  2  kz  ^ pe  j  /p  \  9  ,  2  7  \ 

1  H - —  1  TT^T  +  _  0  ^3)  -QP  +  UpeJo  (fy?) 


^ccJ  9f2 


Id;)2 

k2  2ul  3 


(84) 


c  1  a  \  9 
£(/°  >3  dt 


-  -  r  1  ujce  .  .  dSn 

Eq±  x  VSn  4 - E0  cos  (u>0t.)  — — , 

iz  uJo  oz 


Q 2 

—  +o;4  )  m  =uliJ0{b)ve, 


(85) 


where  for  63  we  have 


63  = 


ul  - 


rnxJn  V  a;.2 


Eo  x  fe 


,i  2 


0  V  “'ce  L 


+  klE% 


2  M)2> 


(86) 


and  J0  is  the  zeroth-order  Bessel  function. 

Assuming  the  following  characteristic  frequency  range 


kVfi  <?C  12  <C  kzVre, 


(87) 


where  Vre  and  Vn  are  the  electron  and  ion  thermal  velocities,  we  derive  an  expression  for 
the  low  frequency  ion  acoustic  oscillation.  Assuming  that  the  electrons  are  magnetized,  the 


52 


low  frequency  electron  oscillations  can  be  described  by  the  following  equation  of  motion, 


4"-^r+mMv«}-  <88> 

where  the  electron  inertia  term  is  omitted  and  <5\1/  is  the  electrostatic  potential  of  the  para¬ 
metrically  excited  ion  acoustic  oscillations,  and  the  angle  brackets  around  the  last  term 
indicate  averaging  over  the  fast  time  scale  cjq  1.  The  Poisson  equation  for  T  is 


A'P  =  47re  ( Sne  —  6rii) , 


(89) 


where  A  is  the  Laplacian  operator. 

Noting  that  the  ion  motion  is  unmagnetized,  we  obtain  the  following  nonlinear  differential 
equation  for  the  density  perturbations  in  the  low  frequency  oscillations, 


+  (Uoihablaj2  -  V?A 


^  =  n™A((V'A)vez), 


(90) 


where  Vs  is  the  sound  speed. 

We  can  now  use  Eqs.  (83),  (84),  and  (88)  to  obtain  the  dispersion  relation  for  the  para¬ 
metric  instability  of  the  oblique  resonance  oscillations.  First,  using  a  standard  procedure,  we 
introduce  two  high  frequency  oscillations  nonlinearly  excited  by  the  parametric  instability 
with  frequencies  u/j  and  uj2  and  a  low  frequency  ion-acoustic  perturbation  at  frequency  il, 


Wei  =  ^lfrexp  i  (kr—u)\t^  ,  nt\  =  Blkexp  i  (kr  —  u  ttj 


(91) 


We2  =  Mk  exp 


^ kr  -f- 


rii 2  =  Bik  exp 


(kr  +  uj2tj 


(92) 
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(93) 


Sn  =  Srik  exp 


Next  we  substitute  Eqs.  (91)-(93)  into  Eqs.  (84),  (85),  and  (90)  to  obtain  the  disper¬ 
sion  relation  for  simultaneous  parametric  excitation  of  two  high  frequency  oscillations  with 
frequencies  given  by 


10, 


wl,2  — 


pi 


1  +<■& 


/K 


1  + 


U). 


pe 


u 


(94) 


and  a  low-frequency  ion-acoustic  frequency  perturbation.  The  dispersion  relation  has  the 
following  form 


1  e2  a ?LH 


U ce  Wl 


^(&) 


k  x  En 


s2  -  a2' 


(95) 


where  (5  =  u>o  —  aq  and  a,  =  e2  E$  /  M2lJq. 

It  is  necessary  to  point  out  that  we  are  interested  in  the  supersonic  regime  of  parametric 
instability  excitation,  i.e.  when 

ft  »  kVs.  (96) 


u;0  —  oj\  <  0,  Aperiodic  Instability 


For  S  <  0,  aperiodic  instability  takes  place  with  a  peak  growth  rate  of 


Ifmax.aper  —  ^ LH 


1  A/  O-’pe  j2(L\  E2 
4  m  uj2e  0  Anmnoc2 


1/3 


(97) 


54 


To  determine  the  plasma  temperature,  we  take  into  account  that  the  effective  electron 
collision  frequency  is  of  the  order  of  the  instability  growth  rate,  vejf  ~  7.  In  the  vicinity  of 
the  heated  plasma,  condition  (92)  holds,  and  electrons  drift  along  the  magnetic  held  lines 


with  sound  speed  Cs.  Here  the  energy  balance  can  be  expressed  as 


c4("t)= 


iO 'p 

i  +  Hr 

^ ce 


Veff 


Eo 


Sir 


(98) 


Substituting  Eq.  (97)  into  Eq.  (98),  the  following  estimate  of  an  effective  plasma  temperature 
increase  in  the  heated  region  can  be  obtained, 


E0 


8/9 


47mm  c2 


(99) 


/ 


where  L  is  the  characteristic  length  of  the  heated  region. 


cuq  —  <jj\  >0,  Periodic  Instability 

For  8  >  0,  the  temperature  enhancement  can  be  found  analogously.  The  dispersion 
relation  now  takes  the  form 


S!2  +  -  eel  = 

2  4  nr  u/f,  U\ 


k  x  E0 


12  8 


82-W 


(109) 


For  8  >  0,  periodic  instability  occurs  with  a  peak  growth  rate 


^n/3 

'Umax, per  —  ^ 


1  Mwi  El 


1/3 
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/ 2(1) 

16  m  uj2e  0  Aumnoc2 


(101) 
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The  peak  growth  rate  of  the  parametric  instability  can  be  written  in  the  form 


{kae)3/2  ulh  , 


(102) 


where  k  ~  u)pe/c  and  ae  is  the  displacement  of  an  electron  during  the  period  of  oscillations 
due  to  the  drift  in  the  wave  field. 

After  the  same  sequence  of  substitutions  as  in  the  case  of  periodic  instability,  the  rise  in 
plasma  temperature  is  found  to  be 


/  1  \  2/3  / 
T  ~  ( -< ^lhL'/Mjjic2  )  l 


2/9 


2 

E0 

\ 


8/9 


An  nine? 


(103) 
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Equations  (99)  and  (103)  are  similar  in  form  and  predict  similar  temperature  increases. 

As  an  example,  we  consider  a  current  running  in  the  antenna  Io  ~  1  A,  the  plasma  density 
inside  the  sheath  of  10"  cm-3,  and  external  magnetic  field  of  1  T.  We  can  estimate  the  field 
excited  in  the  vicinity  of  the  antenna  to  be  Eo  ~  0.G  statvolt/cin.  From  Eq.  (103)  we  can 
estimate  the  increase  in  plasma  temperature  due  to  the  presence  of  the  electric  field  excited 
by  the  antenna  to  be  T  ~  7  eV.  This  rough  estimate  confirms  that  indeed  the  plasma  in 
the  vicinity  of  the  VLF  antenna  can  undergo  a  significant  increase  in  temperature.  This 
in  turn  can  strongly  influence  the  dielectric  properties  of  the  plasma  and  lead  to  enhanced 
disturbance  of  the  received  signals. 


V.  DISCUSSION  AND  SUMMARY 

The  analysis  of  the  conversion  between  EM  and  EAWs  has  heen  carried  out  for  a  sharp 
plasma  boundary,  and  the  results  demonstrate  reasonable  power  for  the  transmitted  EM 
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waves.  Our  next  step  is  to  evaluate  the  feasibility  of  carrying  out  EM  PIC  simulations  to 
model  the  generation  and  propagation  of  the  EAWs  and  their  transformation  to  EM  waves 
at  the  plasma-vacuum  interface. 

In  a  separate  report  [12],  we  expanded  on  our  previous  analysis  of  the  magnetic  win¬ 
dow  analysis,  adding  revised  2D  PIC  simulations  and  new  3D  PIC  simulations  of  the  wave 
propagating  in  a  microwave  horn  and  the  interaction  with  a  plasma  layer. 

Our  initial  analysis  of  the  ReCoinm  scheme  18.  28.  29],  where  applied  electric  and  mag¬ 
netic  fields  are  used  to  open  a  frequency-space  window  for  EM  wave  propagation  through 
the  plasma  sheath,  is  presented  in  Sec.  III.  We  have  established  that  the  ion  distribution  in 
the  cross-field  dipole  (in  the  plasma-optic  regime)  is  altered  significantly  only  on  diffusion 
time  scales,  for  which  r  oc  L7/C,  where  q  ~  1  —  2.  For  plasma  parameters  believed  to  be 
consistent  with  the  radio  blackout  regime,  r  is  of  the  order  of  tens  of  /is,  when  L  on  the 
order  of  a  few  cm,  and  V  is  on  the  order  of  100  V.  For  a  realistic  ReComm  device  allowing 
communication  in  the  L-band  of  the  radio  spectrum,  the  diode  length  L  would  have  to  be 
scaled  up  to  something  on  the  order  of  1  m.  The  magnitude  of  the  applied  voltage  V  is, 
of  course,  limited  by  practical  considerations.  The  utility  of  the  ReComm  scheme  depends 
upon  being  able  to  alter  the  plasma  density  profile  on  a  timescale  which  is  short  compared 
to  the  time  it  takes  the  hypersonically  flowing  plasma  to  pass  by  the  diode.  So  it  question¬ 
able  whether  r  can  be  made  small  enough  to  be  useful  in  alleviating  the  blackout  problem. 
As  mentioned  above,  scaling  the  electrode  separation  up  to  lengths  of  1  m  would  require 
dropping  the  magnetic  field  strength  to  remain  in  the  plasma  optic  regime.  This  will  alter 
the  magnetization  (ratio  of  cyclotron  and  collision  frequencies)  of  the  plasma  species  and 
could  alter  the  scaling  at  longer  diode  lengths. 

We  have  shown  that  an  ambipolar  diffusion  theory  of  the  cross-field  diode  can  describe 
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qualitatively  the  results  of  the  simulation.  However,  plasma  energy  transport  and  charge 
separation  effects,  neglected  in  the  theory,  play  an  important  role  in  determining  the  details 
of  the  plasma  evolution.  We  have  shown  in  our  baseline  simulation  that  ion  temperatures 
can  reach  values  of  several  tens  of  eV  on  long  time  scales  t  ~  r.  This  suggests  that  Ohmic 
heating  of  the  plasma  may  lead  to  further  ionization.  To  examine  such  effects  would  require 
replacing  the  simple  Drude  collision  model  used  in  this  report  with  a  more  realistic  chemistry 
model.  Although  it  is  probably  more  important  to  establish  initially  that  r  can  be  made 
small  enough  for  practical  use. 

Another  concern  is  that  the  direction  in  which  the  plasma  is  shifted  by  the  applied  held 
depends  sensitively  on  the  relative  mobility  of  the  electrons  and  ions.  To  determine  the 
relative  mobility,  it  is  necessary  to  have  reliable  information  about  the  electron  and  ion 
temperatures  and  collision  frequencies  in  the  sheath  surrounding  the  vehicle.  Our  survey  of 
the  literature  suggests  that  the  existence  of  this  information  is  still  incomplete. 

We  note  that  the  ReComm  scheme  is  based  only  on  altering  the  local  plasma  density  in 
the  vicinity  of  the  antenna.  To  do  this  requires  both  an  applied  magnetic  hold  and  a  high 
voltage  diode.  Moreover  the  effect  may  occur  only  on  prohibitively  long  time  scales.  By 
contrast  the  magnetic  window  approach  [30],  which  utilizes  the  magnet  but  does  not  require 
the  diode,  does  not  require  any  slow  bulk  motion  of  the  plasma,  but  instead  alters  the 
transparency  of  the  existing  plasma  to  the  RF  signal.  But  the  effectiveness  of  this  alternate 
scheme  depends  strongly  on  the  plasma  collisionality  [12],  which,  as  noted  above,  is  not  well 
known . 

In  all  communication  scenarios  under  consideration,  the  computational  models  will  be 
extended  to  include,  with  increasing  levels  of  complexity,  realistic  gas  and  plasma  properties 
associated  with  hypersonic  flight  [31  33].  The  goal  of  this  integrated  theoretical  analysis  is 
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the  formulation  of  “operational-windows”  in  which  physical  constraints  associated  with  each 
communication  scheme  can  be  parameterized.  This  analysis  should  significantly  impact  the 
down-select  of  these  schemes  for  eventual  deployment  on  hypersonic  vehicles. 
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APPENDIX  A:  NUMERICAL  CALCULATION  OF  THE  TRANSFORMATION 
COEFFICIENT 


The  following  C  computer  program  was  developed  by  Dr.  Vladimir  Sotnikov  to  numeri¬ 
cally  calculate  the  wave  transformation  coefficient.  Electronic  copies  of  the  program  can  be 
obtained  directly  from  Dr.  Sotnikov  (sotnikov@unr.edu). 

\  /***************************************************/ 


/*  03-16-08  */ 

/*  */ 

/*  EAW_TR  */ 

5  /*  Electron  Acoustic  wave  transformation  */ 

/*  */ 


/*  this  program  creates  the  data  array  output  */ 
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/*  with  2  columns  (theta3,  WT) .  theta3  is  the  */ 
/*  propagation  angle  of  the  outgoing  EM  wave.  */ 
/*  WT  is  the  coefficient  of  transformation  */ 
/*  calculated  as  the  ratio  of  energy  fluxes  */ 
/*  (Pointing  vectors)  of  outgoing  EM  wave  to  */ 
/*  the  incident  EAW.  Description  is  in  reports  */ 
/*  for  January,  February  and  March  2008.  */ 
/*  */ 
/*  */ 


/************+******★****★*★*★★*********************/ 

#include  <iostream.h> 

#include  <stdio.h> 

#include  <math.h> 

main() 

{ 

int  i ; 

float  Tec,  Teh,  nec,  neh; 

float  omega,  omegapec,  omegapeh; 

float  omega2,  omegapec2,  omegapeh2; 

float  c,  pi,  VTec,  VTeh,  lambdas; 

float  absepsc,  absepstot,  G; 

float  rDec,  rDec2,  rDeh,  rDeh2,  RD,  RD2; 

float  al,a2,a3,  invrDc2,  invrDh2,  delta,  ql ; 

float  invrDec2,  invrDeh2,  invrDeh; 

Tec  =  0.2; 

Teh  =  20.0; 
nec  =  0.86ell; 
neh  =  0.05ell; 
pi  =  3.14; 
c  =  3el0; 

VTec  =  4. 19e7*sqrt (Tec) ; 

VTeh  =  4. 19e7*sqrt (Teh) ; 
omegapec  -  5. 64e4*sqrt (nec) ; 
omegapec2  =  omegapec* omegapec; 
omega  =  0 . l*omegapec; 
omega2  =  omega*omega; 
omegapeh  =  5 . 64e4*sqrt (neh) ; 
omegapeh2  =  ome gapeh* omegapeh ; 
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rDec  =  VTec/omegapec; 
rDeh  =  VTeh/omegapeh ; 
rDeh2  =rDeh*rDeh; 
rDec2  =  rDec*rDec; 
invrDeh  =  1.0/rDeh; 
cout  <<  invrDeh; 
invrDec2  =  1.0/rDec2; 
invrDeh2  =  1.0/rDeh2; 

RD2  =  CrDec2*rDeh2)/(rDec2+rDeh2) ; 

//aa  =  sqrt (omegapec2/omega2) ; 

ql  =  invrDeh* (1 .0/sqrt ( (omegapec2/omega2)-l .0)) ; 

//ql  =  invrDeh* (1. 0/sqrt ( (omegapec2/omega2) ) ) ; 
lambdas  =  (2*pi)/ql; 

//lambdas  =  1.0; 

cout  <<  lambdas; 

al  =  VTeh/c ; 

a2  =  omegapec/omegapeh; 

a3  -  omega/ (c*ql) ; 

absepsc  =  (omegapec2/omega2)-l . 0; 

G  =  1 . 0+(RD2/rDeh2)*(l .0/absepsc) ; 

absepstot  -  ( (omegapec2+omegapeh2)/omega2)  -  1.0; 

/*******************************:*************/ 

FILE  *fw ; 

fw  =  f open(Mplasma_parameters . dat" ,  "w"); 

fprintf(fw,  ,,0/of  */0f  70f  \n",  omegapec,  omegapeh,  omega); 
fprintf(fw,  ,,#/0f  7of  ‘/of  ‘/of  \n” ,  VTec,  VTeh,  absepstot,  ql); 

f close (fw) ; 

float  theta3[201]; 
float  WT [201] ; 
float  aa[201]; 
float  bb  [201]  ; 
float  cc[201]; 
float  dd  [201]  ; 

delta  =  (pi/2. 0)/200.0  ; 
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FILE  *fwl; 


fwl  =  fopenO'WT.dat" ,  "w"); 


for  (  i  =  0;  i  <=  200;  i++  ) 

{ 

theta3[i]  =  delta*i; 
aa[i]  = 

16.0*pi*absepsc*al*a2*cos(theta3[i] )*sin(theta3 [i] )*sin(theta3 [i] ) ; 


bb  [i] 


(a3/G) *sin(theta3  [i] ) *sin(theta3 [i] )-absepstot*cos(theta3 [i] ) ; 
dd[i]  =  bb[i]*bb[i]; 

cc[i]  =  absepsc  +  sin(theta3 [i] )*sin(theta3 [i] ) ; 

WT  [i]  =  aa[i]  / (dd [i]  +cc  [i]  )  ; 

//WT [i]  =  ( 16 . 0*pi*absepsc*al*a2*cos(theta3 [i] )*sin(theta3 [i]  )*sin(theta3 [i] ) ) / 

// 

( ( (a3/G) *sin(theta3 [i] )*sin(theta3 [i] )-absepstot*cos(theta3 [i] ) )* 

//  ( (a3/G) *sin(theta3 [i] )*sin(theta3 [i] )-absepstot*cos(theta3 [i] ) ) 


+ 


// 


(absepsc  +  sin(theta3 [i] )*sin(theta3 [i] ) ) )  ; 


fprintf  (fwl ,  "7od  %  f  If  %f  V.f  7d  \n\  i,  theta3[i]  ,  WT[i]); 

> 


f close(fwl) ; 


> 
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APPENDIX  B:  AMBIPOLAR  DIFFUSION  IN  A  CROSSED-FIELD  DIODE 


We  review  the  problem  of  ambipolar  diffusion  in  a  slab,  which  is  treated  in  some  detail 
in  Ref.  [24]  (see  pp  157-163  of  the  reference).  We  expand  the  treatment  to  allow  for  the 
application  of  an  applied  voltage  across  the  electrode  plates  containing  the  plasma. 

The  ambipolar  diffusion  equation  comes  from  considering  continuity  and  momentum 
equations  for  both  ions  and  electrons  in  which  quasi-neutrality  is  assumed  to  hold,  that 
is,  ne  ~  rti.  As  with  the  simulations,  we  assume  that  collisionality  for  each  plasma  species  is 
dominated  by  neutrals,  and  represented  by  a  constant  collision  frequency  vsn ,  where  s  —  i 
or  e,  for  ions  and  electrons,  respectively.  We  also  assume  constant  plasma  temperatures. 
If  the  convective  derivative  in  the  momentum  equations  can  be  neglected,  an  Ohm’s  law 
expression  results 

vs  =  ±/i, E  -  D9V( Inn),  (Bl) 

where  n  is  the  quasi-neutral  plasma  density.  The  top(bottom)  sign  in  Eq.  (Bl)  applies  to 
ions  (electrons).  The  mobility  of  species  s  is  defined  as 


Vs  = 


e 

irtsVsn 


(B2) 


where  we  have  assumed  Z  =  1  for  all  ions,  and  the  diffusion  coefficient  is  given  by 


msvSn ' 


(B3) 
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To  maintain  quasi- neutrality,  it  is  necessary  to  equate  electron  and  ion  fluxes,  n\ *  =  nve. 


From  this  we  obtain  the  ambipolar  electric  field 


E  =  ViV(lnn), 


(B4) 


where  Va  has  dimensions  of  voltage,  and  is  given  by 


Ah  T  /ie 


(B5) 


Note  that  the  sign  of  Va  depends  on  the  relative  magnitudes  of  the  electron  and  ion  diffusion 


coefficient  magnitudes.  When  Va  <  0,  electrons  are  the  more  “mobile”  species  and  the 
ambipolar  held  adjusts  to  slow  down  the  electrons  so  that  their  flux  remains  equal  to  that 
of  the  less  mobile  ions.  When  Va  >  0,  the  ions  are  the  more  mobile  species. 

Inserting  the  common  hux  into  either  the  ion  or  electron  continuity  equation  results  in  a 
diffusion  equation  with  the  ambipolar  diffusion  coefficient 


Ah^e  +  Me-fA 

Ah  +  Me 


(B6) 


Equations  (B‘2)  and  (B3)  were  derived  for  an  unmagnetized  plasma.  As  shown  in  Ref.  [24] 
(pp  169-173),  for  diffusion  perpendicular  to  the  magnetic  held  the  mobility  and  diffusion 
coefficients  for  each  species  are  easily  modified: 


e 


1 


msvsnl  +  (ucs/vsn)2' 


ms^sn  1  +  (Wcs/lAjn)2 


(B7) 
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The  expressions  for  the  ambipolar  voltage  Va  and  diffusion  coefficient  Da  across  Held  lines 
are  still  given  by  Eqs.  (B5)  and  (B6),  as  long  as  y,s  and  Ds  are  calculated  by  Eq.  (B7). 

For  the  parameters  in  Table  I,  c oce/uen  and  I'in/^d  are  both  much  greater  than  one.  In 
this  regime 


K  ~ 


T{  LOpp  Tf  l  A 


Da 


ce 

e  ven 
Te  +  Tj 

T7leU)re 


e  "in 
6  UJci 
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^eri  M cx 


^ °e  IS  m 

Ven  d 
■1 


-1 


(B8) 


Actually,  for  the  parameters  in  Table  I,  ujce/uen  and  Vin/^ci  only  differ  by  about  25%.  In 
this  case,  using  the  initial  plasma  temperatures,  we  can  make  the  simple  approximation 
eVa  ~  —Tej 2  —  5  eV.  Since  Va  <  0,  it  follows  from  Eq.  (B5),  that  De  >  Di,  or  that  the 
electrons  are  the  more  mobile  of  the  two  species. 

We  now  look  for  a  solution  of  the  ambipolar  diffusion  equation  in  one  dimension  (d/dy  — 
d/dz  -  0).  We  consider  the  problem  of  a  plasma  conHned  between  two  conducting  plates  of 
large  extent.  Namely  a  cathode  at  .r  =  0,  and  an  anode  at  x  =  L.  For  electrodes  with  length 
and  width  comparable  to  the  plate  separation  L,  it  is  important  to  consider  diffusion  both 
parallel  and  perpendicular  to  the  magnetic  field.  But  we  neglect  these  multi-dimensional 
effects  for  the  present. 

We  follow  the  treatment  of  Ref.  [24]  and  search  for  a  solution  by  separation  of  variables 


n(x,  t)  =  e  l/TS(x). 


(B9) 
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The  diffusion  equation  is  satisfied  when  S(x )  is  any  linear  combination  of  cos  ( kx )  and 
sin  (kx)  consistent  with  the  boundary  conditions  of  the  problem,  where 

k2  =  1  /(D0t).  (BIO) 

Integrating  Eq.  (B4)  from  the  cathode  surface  at  x  =  0  to  an  arbitrary  field  point  x  yields 
the  electrostatic  potential 


4>(ar)  —  $(0)  =  —Va In  [S(x)/S( 0)]. 


(Bll) 


The  applied  voltage  difference  V  is  given  by 

V  =  4>(0)  -  4>(L). 


(B12) 


We  also  define  the  mid-potential  Vm  as  the  voltage  difference  between  the  center  and  the 
cathode, 

Vm  =  *(L/  2)-<h(0).  (B13) 

If  we  choose  4>(L)  =  0,  Vm  =  $(L/2)  —  V.  From  Eqs.  (Bll),  (B12)  we  hnd 


Sa 

Sk 


=  ev>v\ 


(B14) 


where  Sk  =  5(0)  and  Sa  —  S(L).  This  equation  constitutes  a  boundary  condition  on  the 
solution  for  S(x).  We  can  also  show  [Eqs.  (Bll)  and  (B 13)]  that 


5a- 


~>M 


-  eVM/Va^ 


(B15) 
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where  Sm  =  S(L/ 2).  A  general  solution  for  S(x)  can  then  be  written  in  the  form 


S(x)  =  SK 


,VfVa 


+  ^  cos  [k(x  -  1/2)]  +  (eV/Vo  _  ^  sill  [k(x  -  L/2)] 


cos/cL/2 


sin  kLj  2 


(BIO) 


From  Eqs.  (B14),  (B15),  and  (B16)  we  hnd 


cos  ( kL/2 )  = 


Sk 

2Sm 


(■ ev/K  +  1)  = 


(ev/Va  +  1)  = 


Sk  +  Sa 

2  SM 


(B17) 


This  equation  constitutes  a  general  dispersion  relation  for  the  wavenumber  k. 

We  now  consider  the  special  case  with  no  applied  voltage  V  =  0.  In  this  case  the  density 


values  at  the  cathode  and  anode  are  equal,  Sk  —  Sa,  and  Eq.  (B17)  reduces  to 


cos  (kL/2)  - 


Sk 

Sm 


=  e 


vM/va 


(B18) 


Note  that  a  physical  solution  requires  that  Vm  and  V0  have  the  opposite  sign,  and  that 
Sk  <  Sm,  although  we  expect  the  density  at  the  electrodes  to  be  much  smaller  than  in  the 
middle  of  the  plasma,  Sk  Sm- 

In  Ref.  [24],  the  diffusion  equation  is  solved  with  the  boundary  condition  Sk  =  0.  The 
wavenumber  k  is  then  obtained  from  cos  (kL/2)  =  0.  For  the  lowest  order  mode  kL  =  7r,  or 


r  = 


(B19) 


with  higher  modes  decaying  away  on  time  scales  shorter  than  r.  We  note  that  as  Sk/Sm  — > 
0,  the  mid  potential  Vm  — 1 ►  oo.  This  is  clearly  not  physical,  so  some  small  but  finite  value 
of  Sk/Sm  should  be  retained.  In  this  case,  for  the  lowest  order  mode,  k  <  n/L,  which 
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corresponds  to  a  longer  time  constant  than  predicted  by  Eq.  (B19).  When  Va  <  0  (electrons 
are  held  back  by  the  ambipolar  field),  the  potential  is  concave  down  at  x  =  L/2. 

We  consider  somewhat  more  carefully  here  what  happens  when  Sk  =  0.  In  this  case  the 
density  profile  becomes 

S(x)  =  Sm  cos  [n(x/L  —  1/2)].  (B20) 

We  see  from  Eqs.  (Bl)  and  (B4),  that  this  density  profile  leads  to  infinite  fluid  velocities  at 
the  cathode  and  anode  surfaces.  So  near  the  walls,  the  diffusion  approximation  breaks  down, 
i.e.  it  is  no  longer  appropriate  to  neglect  the  convective  derivative  in  the  fluid  momentum 
equations.  From  Eq.  (B4),  the  ambipolar  electric  field  also  becomes  infinite  at  the  bound¬ 
aries.  However,  the  electric  field  is  derived  by  arguing  that  quasi-neutrality  justifies  equating 
electron  and  ion  fluxes.  But  quasi-neutrality  is  not  maintained  near  the  boundary  surfaces. 
The  scale-length  for  the  non-neutral  region  is  the  plasma  Debye  length,  so  the  treatment 
given  in  Ref.  [24]  gives  a  good  estimate  for  the  density  profile  when  L  is  much  larger  than 
the  Debye  length.  The  actual  physical  value  of  V\i/Va,  however,  depends  on  the  details  of 
the  non-neutral  sheath  physics,  and  cannot  be  directly  calculated  from  a  quasi-neutral  fluid 
theory. 

We  now  consider  the  case  of  a  non-zero  applied  voltage  V .  From  Eq.  (B16),  the  plasma 
center-of-mass  may  be  calculated  to  be 


(x-i/2) 

L 


1  /  ev^a  —  1  \  /  cos  (a)  \ 

2  ye'//v°  +  1  /  \sin2  (a) ) 


(B21) 


where  a  —  kLj 2.  Note  that  in  the  absence  of  an  applied  voltage,  V /Va  — >  0,  the  plasma 
center-of-mass  is  exactly  in  between  the  two  plates,  (x)  =  L/2.  For  a  finite  voltage  V  the 
center-of-mass  may  shift  from  the  center.  Determining  the  magnitude  of  the  shift,  however, 
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requires  solving  the  dispersion  relation  for  a.  From  Eq.  (B17),  this  requires  knowledge 
of  the  two  parameters  V/Va  arid  Vm/Vo.  For  the  baseline  simulation  parameters  given  in 
Table  I,  V  =  100  V,  and  Va  ~  —5  V,  for  which  ev <C  1.  In  this  limit,  the  first  quantity 
in  parenthesis  on  the  right  hand  side  of  Eq.  (B21)  is  approximately  equal  to  -1,  and  the 
dispersion  relation  becomes 

g  ^  a. 

c,os(of)~ — - — .  (B22) 

As  mentioned  above,  the  value  of  the  mid  potential  Vu  cannot  be  determined  by  quasi¬ 
neutral  fluid  theory,  but  depends  on  the  details  of  the  sheath  physics.  Note  that  in  contrast 
to  the  case  with  V  =  0,  a  physical  solution  does  not  require  Vm/Vo  <  0.  In  this  case 
it  sufficient  that  eVM^Va  <  2.  As  Vm/Vo  is  decreased  from  its  theoretical  maximum  of 
+  ln(2),  (x  —  L/2)/L  increases  monotonically  from  —0.25,  to  its  asymptotic  value  of  zero 
for  Vm/Vo  —  1,  while  a  increases  monotonically  from  zero  to  rr/2. 

If  we  estimate  that  Vm/Vo  ~  —1  ( Sk/Sm  ~  0.37),  we  obtain  a  ~  1.39,  and 


(x  ~  L/2) 

L 


~  -0.05, 


(B23) 


that  is,  the  ion  center  of  mass  shifts  by  5%  towards  the  cathode.  Suppose  that  for  different 
values  of  electron  and  ion  temperatures  and  collision  frequencies,  it  occurs  that  Va  >  0. 
but  that  the  magnitude  of  Va  is  still  much  less  than  V,  In  this  case,  the  first  quantity  in 
parenthesis  on  the  right  hand  side  of  Eq.  (B21)  switches  from  ~  —  1  to  ~  +1.  If  we  still 
assume  that  Vm/Vo  ~  —1,  then  the  center-of-mass  shifts  towards  the  anode  side  by  5%, 
instead  of  the  cathode  side.  This  demonstrates  how  the  relative  signs  of  V  and  Va  determine 
the  direction  of  plasma  shift.  In  fact,  we  can  state  the  following  general  rule:  when  Va  is 
negative  (electrons  more  mobile  than  ions),  the  density  distribution  will  skew  towards  the 
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cathode.  When  Va  is  positive  (ions  more  mobile),  the  distribution  skews  towards  the  anode. 
This  is  regardless  of  the  actual  value  of  V\i  (provided  that  it  produces  a  physical  result). 
The  magnitude  of  the  shift  in  center-of-inass  does  of  course  depend  on  V\i. 

This  theory  predicts  that  for  times  t  ~  r,  after  all  higher-order  modes  have  decayed 
away,  ( x )  should  be  independent  of  time  [a  consequence  of  the  ansatz  in  Eq.  (B9)].  This  is 
in  contrast  to  the  simulation  results  which  show  (x)  decreasing  monotonically  with  time,  even 
for  t  >  T.  We  must  recall  that  this  theory  neglects  important  plasma  transport  phenomena. 
For  example,  the  assumption  of  isothermal  plasmas  neglects  Ohmic  heating.  Since  Ohmic 
heating  scales  as  the  square  of  the  electric  field,  we  expect  heating  of  the  plasma  species  to 
be  greater  nearer  the  electrode  where  the  ambipolar  fields  are  largest.  This  will  introduce  a 
temperature  profile  which  is  not  spatially  uniform,  which  violates  one  of  the  key  assumptions 
used  to  derive  the  simple  diffusion  equation.  In  effect,  the  temperature  gradient  will  cause 
variations  in  the  diffusion  coefficient.  Nonetheless,  this  theory  does  give  a  qualitative  (and 
semi-quantitative)  description  of  the  simulation  results  presented  in  this  report. 
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